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ABSTRACT 


Van  Benthem,  Mark  Hilary,  Ph.D.,  Department  of  Chemistry,  College  of  Science  and 
Mathematics,  North  Dakota  State  University,  November  1995.  Chemometric  Analysis  of 
Multimode  Fluorescence  Data  Obtained  With  a  Pulsed  Tunable  Laser.  Adviser:  Dr. 
Gregory  D.  Gillispie. 

This  research  evaluated  the  capabihties  of  various  chemometric  methods  for  analysis  of 
three-mode  fluorescence  data.  Data  were  collected  using  pulsed  Nd:YAG  laser  and  pulsed 
Nd:YAG  laser-pumped  dye  laser  excitation  and  monochromator-PMT-digital  oscilloscope 
detection.  This  apparatus  produced  data  in  the  form  of  nanosecond  scale  time  decay 
profiles  at  numerous  emission  wavelengths  generating  a  wavelength-time  matrix  (WTM). 
Third-order  data  were  produced  by  varying  analyte  concentrations  or  changing  excitation 
wavelength  to  produce  a  time-resolved  excitation-emission  matrix  (TREEM).  Trilinear 
decomposition  (TLD)  and  global  analysis  methods  were  applied  to  a  WTM-concentration 
3-array  and  a  TREEM.  TLD  methods  and  linear  discriminant  analysis  and  classification 
were  performed  on  highly  complicated  data  in  the  form  of  WTMs  of  various  fuels. 

The  three-mode  data  were  decomposed  with  an  eigenanalysis-based  procedure  (EBP); 
three-mode  alternating  least  squares  (3M-ALS),  also  known  as  PARAFAC;  and  three-mode 
nonnegative  alternating  least  squares  (3M-NNALS). 

Various  rank  estimation  procedures  were  evaluated  for  the  two-mode  and  three-mode 
data.  The  efficacy  of  classification  and  clustering  algorithms  applied  to  reduced  forms  of 
three-mode  fuel  fluorescence  data  was  also  demonstrated. 

Analyses  of  the  WTM-concentration  3-array  illustrate  the  capabilities  of  various  rank 
estimation  and  profile  extraction  techniques  with  low-rank  data.  This  four-component 
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WTM-concentration  3-array  was  easily  decomposed  with  TLD  methods.  3M-ALS  and 
3M-NNALS  offered  a  slight  improvement  to  the  EBP  result.  3M-NNALS  offered  a  minor 
computation  time  improvement  over  3M-ALS. 

A  TREEM  of  a  four-component  solution  of  fluorene,  naphthalene,  carbazole,  and 
phenanthrene  in  water  was  measured.  TLD  was  difficult  because  of  extensive  spectral  and 
time-mode  overlap.  The  EBP  and  3M-ALS  were  unable  to  provide  realistic  factor  profiles. 
3M-NNALS  produced  very  good  results,  generating  recognizable  factors  in  one-tenth  of  the 
time  required  by  3M-ALS. 

WTMs  of  fuels  on  soil  matrices  were  analyzed.  TLD  using  an  EBP  was  not  successful, 
yielding  complex  eigenvectors  or  chemically  meaningless  factors.  3M-ALS  performed 
better  than  the  EBP,  but  it  produced  multiple  factor  degeneracies  and  took  a  great  deal  of 
computation  time.  3M-NNALS  performed  the  best  of  the  three  TLD  procedures. 
Computation  time  was  at  least  two  orders  of  magnitude  faster  than  3M-ALS. 
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1.  INTRODUCTION 


Chemical  analysis  provides  information  about  the  composition  of  a  material  (which 
might  be  groundwater  in  a  well,  soil  under  a  filling  station,  an  ingot  of  a  metal  alloy,  etc.) 
to  address  a  concern  someone  has  about  the  material.^  For  example,  is  water  in  the  well 
safe  to  drink?  Is  there  underground  petroleum  contamination  in  the  property  for  sale?  Can 
this  alloy  be  used  in  a  turbine?  Answering  the  fundamental  question  posed  to  the  analyst 
requires  determination  of  one  or  more  chemical  parameters.  The  desired  parameters  are 
sometimes  qualitative  in  nature,  i.e.,  identities  of  chemical  components  (what  heavy  metals 
are  in  the  water,  petroleum  products  in  the  soil,  strategic  metals  in  the  alloy),  but  more 
often  quantitative  information,  i.e.,  the  amounts  of  chemical  components  (how  much  lead  is 
in  the  groundwater,  gasoline  in  the  soil,  or  chromium  in  the  alloy,  etc.)  is  desired.  The  job 
of  the  analytical  chemist  is  to  obtain  the  best  information  in  the  necessary  form  to  answer 
these  questions.  It  is  in  the  interest  of  both  the  analytical  chemist  and  the  client  to  obtain 
accurate  data  in  a  fast,  inexpensive,  and  minimally  invasive  fashion. 

Traditional  analysis  methods  require  collection  of  representative  samples  from  the  bulk 
material  for  subsequent  laboratory  analysis.  Laboratory  analysis  methods  provide  great 
specificity  at  the  expense  of  slow  speed  and  the  requirement  of  extensive  sample 
preparation  and  tedious  data  analysis.  In  addition,  the  precision  and  the  accuracy  of  the 
data  are  often  limited  by  the  quality  of  the  sampling,^  which  may  not  have  been  performed 
by  the  chemist;^  for  environmental  projects,  the  chemist,  in  fact,  rarely  has  any  connection 
with  the  sampling.  Analytical  techniques  that  allow  measurement  of  bulk  materials  in  real 
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time  and  in  situ  would  be  a  windfall  to  analytical  chemistry.  Fiber-optic  spectroscopy 
provides  such  utility  in  a  number  of  situations. 

Fiber-optic  spectroscopy  can  be  extremely  valuable  in  situations  where  sampling  is  very 
difficult,  dangerous,  or  time  critical.  It  may  also  be  indicated  when  analysis  involves 
obtaining  information  about  a  remote  or  otherwise  inaccessible  location,  for  example, 
spectroscopic  monitoring  of  jet  engine  exhaust  in  flight.  One  may  also  use  it  to  monitor 
multiple  stages  of  an  industrial  process  simultaneously.  Safety  while  performing  analysis  is 
another  issue  that  may  dictate  the  use  of  fiber  optics.  Spectroscopic  analysis  of  explosives 
or  radioactive  materials  are  examples  of  analyses  that  may  be  more  safely  performed  at  a 
distance  using  fiber  optics. 

Near-infrared  (NIR)  absorbance  spectroscopy  performed  over  fiber  optics  has  been 
practiced  for  many  years  ^  and  has  reached  a  relatively  advanced  stage  of  commercial 
maturity.  NIR  instruments  are  used  for  industrial  process  monitoring  and  medical 
diagnostics  and  in  environmental  chemistry.  In  applications  for  which  quantitative 
information  is  required  for  species  at  high  concentrations,  it  can  be  a  superior  analysis 
method,  saving  time  and  money.  Because  samples  are  not  taken,  however,  special 
consideration  must  be  given  to  interferences,  i.e.,  absorbance  by  species  other  than  the 
target  analytes. 

Traditional  analytical  methods  provide  opportunities  to  eliminate  sample  interferences 
before  the  actual  data  collection.  Great  effort  is  expended  to  eliminate  interferences  since 
they  can  be  a  dominant  source  of  error.  If  the  analyte  response  caimot  be  separated  from 
that  of  the  interference,  one  must  model  the  influence  of  the  interferences  on  the  analyte 
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response.  Fiber-optic  spectroscopy  and  other  in  situ  techniques,  which  generally  do  not 
allow  sample  manipulation,  are  especially  susceptible  to  interferences  and  matrix  effects. 
The  burden  of  detecting  and  treating  the  interferences  falls  on  the  quality  of  the 
spectroscopic  d^ta  and^  ultimately,  on  the  data  analysis. 

To  correctly  predict  the  concentration  of  a  spectroscopically  active  species,  all 
spectroscopically  active  species  present  in  the  unknown  must  be  included  in  the  calibration 
model.  If  there  are  interfering  species  whose  contributions  were  not  accounted  for  in  the 
calibration  procedure,  the  concentration  determination  will  be  biased.^  For  chemically 

complex  systems,  as  are  often  encountered  in  environmental  analyses,  interferences  can  be 
a  difficult  problem  to  overcome. 

Multivariate  analysis  of  NIR  spectroscopic  data  represents  a  highly  refined  approach  to 
modeling  absorbance  data  and  is  the  source  of  the  technique's  popularity  and  utility. 
However,  because  NIR  (and  other  absorbance  techniques)  are  limited  to  a  single  data 
mode,  viz.,  wavelength,  their  ability  to  predict  unknown  parameters  in  the  presence  of 
mterferences  not  included  in  the  calibration  model  is  limited.  Multimode  data,  of  order^ 
greater  than  that  available  in  absorbance,  can  allow  calibration  and  analysis  even  if  there 
are  interferences  not  modeled  in  the  calibration.^’’ 

Fluorescence  spectroscopy  is  readily  adapted  to  multimode  techniques.®  Moreover,  it  is 
several  orders  of  magnitude  more  sensitive  than  absorbance’  and,  unlike  absorbance 
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to  many  worfe,  the  term  dimension  is  used  in  place  of  the  term  mode,  which  wiU  be  used  here.  The  term  mode 
indicates  the  tonchonality  employed  to  acqm’re  data.  For  example,  wavelength  is  the  mode  utilized  in  absorbance 
sp^troscopy  while  charge  to  mass  ratio  is  the  mode  of  mass  spectral  data,  time  is  the  mode  of  gas  chromatographic 
data,  etc.  Mode  is  explained  in  greater  detail  in  Appendix  B.  ^  ^ 


T^e  meaning  of  the  term  order  is  discussed  in  greater  detail  in  Chapter  3;  however,  consider  it 
of  modes  used  to  collect  data. 


to  represent  the  number 
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techniques,  is  easily  applied  to  highly  scattering  matrices  like  soils.  By  measuring  the 
fluorescence  emission  spectrum  at  different  excitation  wavelengths  (or  equivalently, 
measuring  the  excitation  spectrum  at  different  emission  wavelengths),  one  can  generate  a 
matrix  of  data  commonly  referred  to  as  an  excitation-emission  matrix  (EEM).  The  two 
modes  of  an  EEM  are  the  excitation  wavelength  and  the  emission  wavelength.  Similarly, 
one  can  measure  the  fluorescence  emission  at  various  emission  or  excitation  wavelengths 
as  a  function  of  time  following  pulsed  excitation  to  generate  a  wavelength-time  matrix 
(WTM).  A  series  of  excitation-emission  matrices  obtained  at  various  time  delays  following 
pulsed  excitation  represents  three  data  modes.  Given  the  higher  order  techniques 
potentially  available  and  the  overall  versatility  of  the  technique,  fluorescence  spectroscopy 
promises  to  be  a  superior  spectroscopic  analysis  method  in  many  circumstances. 

Unfortunately,  fluorescence  spectroscopy  has  an  inherently  lower  signal-to-noise  ratio 
(S/N)  than  absorbance  speetroscopy.^  The  ramifications  here  are  lower  precision 
measurements,  hence,  lower  precision  parameter  determination  (e.g.,  concentration).  Of 
course,  in  situations  that  would  otherwise  lead  to  large  sampling  uncertainty  or  great  risk 
associated  with  sampling,  the  precision  trade-offs  of  a  fiber-optic  fluorescence 
measurement  may  be  acceptable.  Again,  one  must  consider  the  type  of  information  desired 
of  the  measurement  (the  data  quality  objectives). 

If  the  analysis  goal  is  quantitative  information,  the  analysis  requires  a  calibration  step 
followed  by  an  unknown  estimation  step.  The  presumption  in  this  case  is  that  one  knows 
the  identity  of  the  species  present  and  desires  the  concentrations.  However,  when  the  goal 
of  the  analysis  is  qualitative  information,  one  is  faced  with  a  more  thorny  problem; 
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especially  if  the  material  has  many  unknown  components  with  similar  spectra  and  there  are 
multiple  matrix  effects,  as  in  environmental  analysis.  The  parameters  of  interest  now 
become  the  relative  proportions  of  the  spectral  intensities  as  a  function  of  wavelength, 
rather  than  their  absolute  magnitudes,  i.e.,  the  spectra  of  the  contaminants.  The  data  must 
be  decomposed  to  yield  the  individual  spectra  of  the  components. 

The  presence  of  noise  in  the  spectral  data  can  make  decomposition  difficult,  leading  to 
ambiguous  results.  Here,  too,  the  multimode  data  provided  by  fluorescence  compensate 
greatly  in  spite  of  the  lower  S/N  compared  to  absorbance  spectroscopy.  Once  again,  the 
analyst  is  faced  with  trade-offs  between  S/N  and  sensitivity  and  higher  data  order.  In  fact, 
three-mode  (and  higher  mode)  data  with  linear  independence  between  the  modes  (so-called 
trilinear  data)  produce  a  unique  decomposition  if  some  other  mild  assumptions  are  met.^^’^^ 
Simulations  of  decomposition  of  trilinear  data  with  various  levels  of  noise  demonstrate  that 
analysis  of  fluorescence  spectra  is  very  encouraging. 

Given  the  many  potential  advantages  of  fluorescence  spectroscopy,  its  increasing  use  in 
the  chemical  analysis  of  complex  systems  is  not  surprising.  Environmental  analysis  has 
been  impacted  immensely  by  the  marriage  of  fiber  optics  and  fluorescence.  The  EPA's  list 
of  the  20  most  significant  hazardous  substances  includes  several  fluorescent  species.*'^ 

Most  common  fuels  contain  significant  levels  of  fluorescent  compounds,  specifically 
benzene  and  its  single  aromatic  ring  derivatives,  naphthalene  and  substituted  naphthalenes, 
and  other  polycyclic  aromatic  hydrocarbons  (PAHs).*^‘^^  Fuels  also  contain  high  levels  of 
aliphatic  hydrocarbons,  which  are  not  fluorescent  (but  could  be  determined  by  NIR). 
However,  in  terms  of  toxicity  and  carcinogenicity,  the  aliphatic  compounds  do  not  pose 
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nearly  the  threat  that  the  aromatic  species  do.**  It  is  fortunate,  then,  that  fluorescence 
selectively  detects  the  species  of  greatest  concern. 

The  collection  of  truly  representative  environmental  samples  from  a  suspected 
contaminated  site  can  be  extremely  challenging  and  expensive}^  particularly  in  the  case  of 
subsurface  samples.  In  the  standard  approach,  the  areas  of  interest  must  be  drilled  or  core 
sampled.  The  installation  cost  of  a  single  monitoring  well  can  exceed  $10,000,  and  this 
money  is  wasted  if  the  well  is  not  placed  in  an  appropriate  location. 

Once  a  well  has  been  installed,  the  samples  must  be  collected,  packaged,  shipped, 
stored,  and  analyzed  judiciously  so  they  are  not  corrupted.  Loss  of  sample  integrity  could 
occur  at  any  or  several  of  these  steps.  Fuel-contaminated  soil  and  water  samples  that 
contain  easily  lost  volatile  and  semi- volatile  components  must  be  handled  with  great  care. 
The  elapsed  time  between  sampling  and  availability  of  the  chemical  analysis  report  may  be 
several  weeks.  To  make  matters  worse,  after  the  drilling  or  core  sampling  is  finished,  a 
substantial  conduit  remains  for  transport  of  the  contaminants  to  other  soil  strata.  Finally, 
analysis  methods  proscribed  by  EPA  to  extract  contaminants  from  soil  involve  the  use  of 
solvents  which  themselves  are  hazardous  materials.  Methods  that  would  allow  results  in 
real  time  with  lower  costs  are  certainly  desirable.  A  promising  method  to  measure  fuel 
contaminants  in  soil  combines  fiber-optic  fluorescence  with  cone  penetrometry. 

Cone  penetrometry  is  a  method  of  determining  geophysical  parameters  of  subsurface 
soil.  It  was  developed  primarily  to  characterize  geology  of  areas  during  oil,  gas,  and 
mineral  exploration.  It  employs  a  hydraulic  press  mounted  in  a  very  heavy  (=20  tons)  truck 
to  drive  sections  of  stainless  steel  pipe  (=5  cm  diameter)  into  the  earth.  The  pipe  contains 
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various  sensors  to  measure  geotechnical  parameters  such  as  sleeve  friction  and  tip 
resistance;  these  parameters  form  the  basis  of  a  soil-type  characterization  scheme.  A 
sapphire  window  in  the  side  of  the  pipe  close  to  the  tip  of  the  probe  allows  optical  access  to 
the  soU.  Fiber  optics  leading  to  and  from  the  window  allow  direct  fluorescence 
examination  of  the  soil.  A  new  device  has  been  developed  specifically  for  this  type  of 
spectroscopy.  The  Rapid  Optical  Screening  Tool  (ROST™,  Dakota  Technologies,  Inc., 
Fargo,  ND)  allows  collection  of  a  fluorescence  time-resolved  excitation-emission  matrix 
(TREEM)  of  petroleum-contaminated  soil  using  fiber-optic  transmission  of  excitation  and 
emission. 

The  ROST™  utilizes  a  frequency  doubled,  Nd:  YAG  laser  pumped-dye  laser  to  generate 
ultraviolet  laser  excitation  which  can  be  passed  into  an  optical  fiber.  The  detector 
incorporates  a  scanning  monochromator,  photomultiplier  tube  (PMT),  and  digital  storage 
oscilloscope  (DSO).  The  output  of  the  DSO  is  channeled  into  a  personal  computer  (PC) 
and  stored  for  processing.  The  coupling  of  ROST™  with  the  cone  penetrometer  makes  a 
powerful  tool  that  can  quickly  produce  multimode  spectra  of  contaminated  soil  to  rapidly 
characterize  a  contaminated  site.  Compared  to  traditional  methods  of  chemical  analysis, 
cone  penetrometry  with  fiber-optic  detection  is  more  cost  effective  and  may  be  more 
accurate. 

The  rationale  for  pursuing  cone  penetrometry-based  fluorescence  analysis  is 
compelling.  Its  strength,  in  addition  to  the  factors  mentioned  above,  lies  in  the  wealth  of 
information  obtained  during  an  analysis.  In  addition,  after  performing  cone  penetrometry. 


the  hole  produced  can  be  grouted  as  the  probe  is  removed.  This  minimizes  the  spread  of 
any  hazardous  material,  unlike  monitoring  wells. 

Rarely  do  sites  contain  a  single  compound  (e.g.,  benzene).  The  fuels  and  residue  from 
other  energy-producing  activities  (e.g.,  manufactured  gas  plants)  represent  very  complex 
mixtures  of  fluorescing  species.  Consequently,  the  spectra  of  the  contaminated  soils  are 
usually  quite  complex.  The  higher  order  inherent  in  fluorescence  data  may  allow  more 
reliable  and  complete  chemical  analysis  than  do  single-mode  techniques. 

Data  obtained  from  cone  penetrometry-based  fluorescence  is  extensive.  Consider 
measuring  fluorescence  intensity  as  a  function  of  excitation  wavelength,  emission 
wavelength,  fluorescence  decay  time,  and  probe  depth  at  various  locations  on  a  site.  It  is 
easy  to  see  how  large  the  data  sets  can  be.  It  is  not  practical  to  perform  analysis,  other  than 
gross  generalizations,  on  such  enormous  data  sets  without  the  use  of  a  computer. 
Fortunately,  there  is  a  relatively  young  branch  of  chemistry  designed  to  deal  with  data 
analysis  of  large  and  complicated  data  sets:  chemometrics. 

Chemometrics  has  been  defined  as  using  chemical  principles  and  mathematical  and 
statistical  methods  to  interpret  and  predict  chemical  data.  Much  of  chemometrics  is  an 
offshoot  of  developments  in  applied  statistics  used  by  psychometricians  and 
econometricians  to  analyze  data  related  to  human  behavior  and  the  performance  of  markets, 
respectively.^^  The  general  idea  is  to  take  a  set  of  data  pertaining  to  a  group  of  subjects  and 
draw  inferences  about  the  nature  of  the  group.  For  example,  a  psychologist  might  want  to 
analyze  infant  behavior  in  relation  to  different  stimuli  to  leam  about  language  development. 
It  is  also  valuable  to  process  the  data  to  find  a  predictive  method  for  subjects  within  or 
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outside  of  the  group.  Economists,  and  almost  everyone  else,  would  like  to  be  able  to 
accurately  predict  when  a  stock  price  wUl  go  up  and  when  it  will  go  down.  It  may  be 
apparent  to  the  reader  that  some  humans  have  exceptional  methods  of  “data  analysis” 
without  the  use  of  computers  and  complex  algorithms.  Given  the  proper  training  and/or 
experience,  people  can  be  veiy  effective  at  detecting  trends,  recognizing  patterns, 
estimating  proportions,  and  predicting  events.  Certainly,  there  are  very  astute  m^ket 
speculators  who  make  a  very  good  living  simply  using  their  “feel”  about  the  market,  just  as 
a  skilled  spectroscopist  can  identify  a  compound  in  a  mixture  and  accurately  estimate  its 
concentration  by  briefly  examining  a  spectrum. 

Analytical  chemistry,  though,  demands  much  more  of  analysis  than  a  rough  guess,  no 
matter  how  accurate.  Chemometrics  provides  the  analytical  chemist  with  numerical  tools 
capable  of  dealing  with  large  and  complicated  sets  of  data  and  reducing  them  to  a  useful 
and  desired  form. 

Data  reduction  is  the  hallmark  of  a  good  chemometric  analysis  routine.  From  large  sets 
of  fluorescence  data,  one  seeks  the  nature  and  quantity  of  the  species  giving  rise  to  the 
fluorescence.  These  may  be  from  contaminants  in  soil  and  groundwater  or  some  other 
source.  The  objective  here,  and  the  subject  of  this  thesis,  was  to  examine  chemometric 
methods  that  could  assist  an  analyst  in  determining  the  identities  of  luminescent  species  in 
solutions  and  mixtures  or  the  identity  of  a  composite  mixture  (e.g.,  a  fuel)  using 
fluorescence  data  obtained  from  the  ROST™  or  a  similar  instrument. 

Chapter  2  of  this  work  reviews  various  properties  of  fluorescence  with  particular 
attention  paid  to  those  aspects  of  fluorescence  exploited  by  the  ROST™.  Chapter  3 
addresses  the  chemometric  methods  used  to  decompose  ROST™  acquired  data.  The 
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principles  of  factor  analysis  are  discussed  to  identify  or  quantify  unknowns.  Methods  of 
separating  collections  of  data  into  groups  using  pattern  recognition  and  classification 
schemes  are  the  subject  of  Chapter  4.  Results  of  the  research  conducted  and 
complementary  discussion  are  presented  in  Chapter  5.  Chapter  6  incorporates  conclusions 
about  this  research  and  recommendations  for  future  work.  For  those  not  familiar  with 
multivariate  analysis  and  matrices,  Appendix  A  compiles  a  selection  of  definitions  and 
explanations  of  concepts  covered  throughout  the  thesis.  Appendixes  B  and  C  also  contain 
background  material  relevant  to  the  subject  of  higher  order  factor  analysis. 
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2.  FLUORESCENCE  REVIEW 


2.1  Fluorescence  Principles 

Fluorescence  in  organic  molecules  is  often  related  to  the  presence  of  delocalized  k 
electrons, most  notably  in  aromatic  hydrocarbons  and  their  derivatives.  In  the  discussion 
of  fluorescence  principles,  we  shall  use  the  language  and  nomenclature  conventions 
relevant  to  the  fluorescence  of  aromatic  molecules. 

Absorption  of  a  photon  by  a  molecule  in  the  ground  electronic  state.  So,  promotes  the 
molecule  into  one  of  the  states,  of  the  excited  electronic  state  manifold.  The  integrated 
strength  of  the  electronic  transition  is  governed  by  an  integral  involving  the  electronic 
wavefunctions,  while  the  shape  of  the  excitation  spectrum  is  determined  by  the 
Franck-Condon  factors.  In  general,  electronically  excited  molecules  in  the  condensed 
phase  rapidly  relax  (~10’^^  s)  to  a  Boltzmann  distribution  of  vibrational  levels  in  the  first 
excited  singlet  state,  Si  via  the  processes  of  internal  conversion  (Sj  Sj )  and  vibrational 
relaxation  ( v' »  0  -»  v'  ~  0).  The  absence  of  appreciable  fluorescence  from  states  higher 
in  energy  than  Si  is  known  as  Kasha’s  rule.^'* 

Deactivation  from  the  Si  state  occurs  through  a  variety  of  meehanisms.  The  radiative 
Si  to  So  deactivation  process  is  referred  to  as  fluorescence.  Others  processes  include 
internal  conversion  (S,  — >  Sq  ),  external  conversion  (energy  transfer  and  quenching),  and 
intersystem  crossing  (singlet  to  triplet  conversion).  The  role  of  the  Franck-Condon  factors 
in  determining  the  shape  of  the  emission  spectrum  is  similar  to  the  way  that  they  affect  the 
excitation  spectrum.  After  returning  to  So,  the  molecule  vibrationally  relaxes  to  a 
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Boltzmann  distribution  of  vibrational  levels.  These  collective  phenomena  have  two 
important  consequences:  (1)  for  a  given  molecule,  the  shape  of  the  fluorescence  emission 
spectrum  (i.e.,  the  relative  intensity  distribution  as  a  function  of  wavelength)  is  independent 
of  the  excitation  wavelength;  and  (2)  the  shape  of  the  fluorescence  excitation  spectrum  is 
independent  of  the  emission  wavelength  monitored  during  its  acquisition. 

The  depopulation  rate  of  the  excited  state  also  has  important  consequences.  The 
differential  equation  for  the  disappearance  of  Si  after  photoexcitation  is 

-%^  =  *,[S,]+t„[S,]+t,[S,]+t,[S,lQ],  (2.1) 

where  [S,  ]  is  the  concentration  of  molecules  in  the  first  excited  state,  is  the  first-order 
rate  constant  for  radiative  transition,  is  the  first-order  rate  constant  for  non-radiative 
transitions  combining  internal  conversion  and  intersystem  crossing,  kp  is  the  first-order 
rate  constant  for  unimolecular  photochemistry,  is  the  second-order  rate  constant  for 
quenching,  and  [Q]  is  the  concentration  of  a  quenching  species.  The  rate  constants  in  (2.1) 
are  usually  combined  into  a  pseudo-first-order  rate  constant,  viz., 

k  =  k,+k„^+kp+  kg  [Q] .  (2.2) 

Thus,  the  fluorescence  follows  first-order  kinetics  and  its  intensity  decays  exponentially  in 
time.  The  fluorescence  lifetime,  T ,  is  the  time  required  for  the  fluorescence  intensity  to  fall 
to  1/e  of  its  initial  value  following  pulsed  excitation  of  short  duration.  The  relationship  of 
lifetime  to  the  rate  constant  is  a  reciprocal  one,  i.e., 

x  =  X-  (2-3) 

Fluorescence  lifetimes  are  generally  in  the  range  of  10'^  s  - 10'^  s.^  Lifetimes  are  affected 
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by  a  variety  of  factors,  including  solvent,  other  solutes,  and  dissolved  oxygen.  There  are 
two  standard  methods  in  use  today  for  measuring  fluorescence  lifetimes:  pulsed  or 
time-resolved  techniques  and  harmonic  or  phase-modulated  techniques.^^’^^ 

2.2  Lifetime  Measurement  Techniques 

The  nature  of  excited  state  lifetime  determinations  and  the  pitfalls  associated  with  the 
different  techniques  are  discussed  in  numerous  sources,  including  a  monograph  by 
Demas.^^  The  two  principal  approaches,  pulsed  and  phase-modulated  techniques,  are 
discussed  in  this  section,  with  emphasis  on  the  former  since  it  was  the  one  used  in  this 
research. 

It  is  convenient  to  relate  lifetime  determination  to  signal  processing  theory  since 
fluorescence  decay  is,  in  reality,  a  time  domain  signal.  Consider  a  fluorescing  species  as  an 
analog  system,  like  an  electronic  black  box.  The  excitation  photons  can  be  viewed  as 
simply  an  input,  e(t),  to  the  system.  The  output  of  the  system  is  the  fluorescence. y(r),  i.e., 

fit)  =  e(t)*y(t),  (2.4) 

where  *  denotes  the  mathematical  operation  of  convolution,  and  y(t)  is  the  impulse 
response  of  the  fluorescence.  Because  these  time-dependent  signals  are  measured  with  a 
detector  and  electronic  signal  processors  that  impose  their  own  time  dependence,  additional 
terms  must  be  incorporated  into  the  equation.  The  system  now  consists  of  the  solution 
containing  the  species  of  interest:  the  detector,  whose  impulse  response  is  d(ty,  and 
processing  electronics,  whose  impulse  response  is  l(t).  The  instrumentally  observed 
fluorescence  is 

/  (0  =  e(t)*  y(t)*  .  (2.5) 
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Since  the  fundamentally  interesting  term  here  is  the  impulse  response  of  the  species 
decay,  i.e.,  y(f),  terms  in  (2.5)  can  be  combined  such  that 

=  (2.6) 

with 

E(t)  =  e(t)*d(t)*l(t).  (2.7) 

E(t)  is  called  the  instrument  response  function,  and  it  can  be  defined  as  the  input  to  the 
system,  y(t) ,  whose  output  is  f(t) .  The  convolution  in  (2.5)  is  defined  by 

/  (0  =  E  ix)y(t  —  x)dx  =  E{t-  x)y{x)dx .  (2.8) 

Since  the  fluorescence  derives  from  an  actual  physical  process,  the  system  is  causal,  i.e.,  the 
signal  at  any  time  does  not  depend  on  future  values  of  the  input  to  the  system.^^  In  fact, 
since  the  impulse  response  is  a  physical  (therefore,  causal)  system  and  y{t)  is  a  causal 
system,  the  result  is 

/(O  =yE{x)y{t  -  x)dx ,  (2.9) 

for  all  t>0  and  zero  for  t<0.  Since  y(t)  has  the  form 

y(0  =  e-*',  (2.10) 

for  all  f  >  0  and  zero  for  t  <  0 ,  (2.8)  becomes 

fit)  =  e-^j‘^E(x)e^dx,  (2.11) 

for  alk  >  0  and  zero  for  t<0,  which  is  the  common  form  for  the  convolution  of  an 
excitation  with  an  exponential  decay. 
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It  is  important  to  note  that  E(t)  is  not  determined  by  just  the  excitation  pulse  shape 
itself.  Rather,  it  reflects  the  response  of  the  entire  system  except  for  the  response  of  the 
fluorophore  itself.  The  formulation  of  the  relationship  in  (2.10)  is  important  when 
considering  how  to  extract  the  fluorescence  lifetime  from  the  data  provided  by  pulsed  and 
phase-modulated  techniques. 

2.2.1  Pulsed  Techniques 

In  the  pulsed  methods,  one  attempts  to  measure  the  lifetime  directly  by  exciting  a 
solution  containing  the  species  of  interest  with  a  short  duration  burst  of  photons.  The 
emission  intensity  is  observed  as  a  function  of  time  following  the  excitation.  Ideally,  the 
response  time  of  the  detector  and  signal  processor  are  much  shorter  than  the  decay,  in 
which  case  one  observes  the  true  decay.  Complications  of  this  simple  process  occur 
because  of  three  factors:  finite  excitation  pulse  width,  i.e.,  e(t);  the  detection  response  time, 
i.e.,  d(t)  with  l(t);  and  experimental  noise. 

There  are  a  variety  of  excitation  sources  one  can  use  to  generate  intense  short  pulses  of 
UV  and  visible  light.  Three  of  the  more  popular  sources  are  (1)  gas-filled  discharge  lamps, 
which  operate  at  a  few  kilohertz  (KHz)  with  pulse  durations  of  a  few  nanoseconds;  (2)  low- 
repetition  rate  (-100  Hz)  excimer,  nitrogen,  and  Nd:YAG  lasers  with  typical  pulse 
durations  of  a  few  nanoseconds;  and  (3)  mode-locked  and  cavity  dtunped  lasers  with  high- 
repetition  rates  (up  to  80  MHz)  and  pulse  durations  in  the  picosecond  range.^® 

Of  course,  the  selection  of  excitation  source  may  depend  upon  the  fluorophore  under 
investigation,  but  it  is  also  affected  by  the  choice  of  detection  scheme.  Ideally,  the  pulse 
duration  is  chosen  to  be  much  shorter  than  the  lifetime,  in  which  case  one  could  consider 
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e{t)  a  Dirac  delta  function,  8(t).  If  this  is  so,  then  one  need  only  consider  the  response  of 
the  detection  system  in  extracting  the  lifetime.  However,  if  the  pulse  duration  is  not  a  delta 
function,  the  data  will  be  smeared,  and  deconvolution  is  necessary. 

Pulsed  detection  schemes  commonly  used  in  measuring  lifetimes  in  the  nanosecond  and 
subnanosecond  regime  include  (1)  time-correlated  single-photon  counting  (TCSPC)  which 
can  measure  picosecond  lifetimes,  (2)  streak  cameras  which  can  also  measure  picosecond 
lifetimes,  and  (3)  oscilloscopes. 

TCSPC  employs  two  PMTs,  signal  amplifiers,  constant  fraction  discriminators  set  up  in 
parallel,  a  time-to-amplitude  converter  (TAC),  and  a  multichannel  analyzer  (MCA).^®  The 
excitation  pulse  simultaneously  illuminates  the  sample  and  the  triggering  PMT.  The 
triggering  PMT  detects  the  excitation  pulse  and  sends  an  amplified  “start  signal”  to  the 
TAC  to  start  the  ramp  generator  of  the  TAC.  If  a  fluorescence  event  is  detected  by  the 
second  PMT,  a  stop  pulse  is  generated,  and  the  TAC  converts  the  ramp  value  to  a  time. 

The  TAC  is  programmed  to  stop  at  a  predetermined  value  if  no  event  occurs.  If  an  event 
occurs,  the  MCA  places  a  count  into  the  appropriate  memory  channel,  thus  building  a 
histogram  of  events  as  a  function  of  time.  If  multiple  photon  events  reach  the  second  PMT, 
the  TAC  will  only  see  the  first  one  and  miss  subsequent  ones.  To  prevent  this  so-called 
“photon  pileup,”  an  average  count  rate  of  0.01  to  0.05  per  pulse  is  desirable.^” 

TCSPC  has  many  advantages.  It  is  one  of  the  most  sensitive  lifetime  approaches  owing 
to  its  need  for  low  intensity  fluorescence.  It  can  measure  decays  much  shorter  than  its 
system  response  would  indicate.  This  results  from  the  triggering  and  counting 
mechanisms.  The  start  and  stop  signals  are  taken  from  the  leading  edge  of  the  excitation 
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pulse  and  fluorescence,  respectively,  as  detected  by  the  PMT  and  determined  by  the  CFD. 
Thus,  a  PMT  with  a  several  nanosecond  response  time  may  be  used  to  measure  a 
nanosecond  lifetime  directly  (provided,  of  course,  the  excitation  source  is  a  diff).  Since 
TCSPC  is  a  counting  experiment,  it  follows  well-defined  Poisson  statistics,  so  error 
analysis  is  greatly  simplified?^  TCSPC  also  has  a  large  linear  dynamic  range,  i.e.,  it  has  the 
ability  to  measure  a  wide  (time)  range  of  decay  processes. 

The  major  disadvantage  in  TCSPC  is  the  time  required  to  collect  a  profile.  For  good 
statistics,  the  recommended  minimum  of  counts  in  the  peak  channel  is  100,000.^°  With 
low  repetition  rate  sources,  acquisition  of  a  decay  curve  can  take  many  minutes  or  even 
hours.  Mode-locked  lasers  are  favored  for  TCSPC  owing  to  their  high  repetition  rates. 

Streak  cameras  can  record  very  fast  processes  in  real  time  by  converting  time  domain 
information  into  spatial  information.  The  process  is  as  follows:  (1)  Incident  photons  are 
converted  by  the  photoelectric  effect  to  electrons  at  a  photocathode.  (2)  The  photoelectrons 
are  accelerated  by  a  potential  difference  into  an  electron  deflection  field,  which  is  varied 
linearly  in  time.  (3)  The  photoelectron  beam  is  swept  across  a  microchannel  plate,  and 
individual  electrons  are  amplified  within  the  microchannels  before  striking  a  phosphor 
screen.  (4)  The  image  from  the  phosphor  is  recorded  by  a  vidicon  or  CCD.  Streak  cameras 
provide  very  high  temporal  resolution,  having  been  extended  into  the  femtosecond  range.^^ 
However,  it  is  difficult  to  get  both  high  temporal  resolution  (<10  ps)  and  long  observation 
times  (>5  ns)  simultaneously.  Streak  cameras  are  also  very  expensive. 

Oscilloscopes  are,  in  general,  the  most  inexpensive  of  the  methods  discussed  here  to 
measure  lifetimes.  In  most  oscilloscope-based  lifetime  studies,  the  anode  signal  from  a 
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PMT  is  directly  connected  to  an  input  channel  of  the  scope.  Oscilloscopes  come  in  two 
basic  varieties,  analog  and  digital.  Analog  oscilloscopes  are  rarely  used  anymore  for  short 
lifetime  (nanosecond  regime)  measurements,  owing  to  the  tediousness  of  the  data 
processing.  Formerly,  images  of  the  analog  oscilloscope  were  captured  photographically 
and  the  decay  curves  digitized  manually  for  subsequent  mathematical  processing.  Digital 
oscilloscopes  or  digital  storage  oscilloscopes  (DSO)  represent  the  state  of  the  art  in 
oscilloscope  technology,  and  many  are  capable  of  recording  transient  signals  in  real-time 
or  averaging  repetitive  transients. 

DSOs  acquire  analog  signals,  the  PMT  output,  by  sampling  and  digitizing.  As  the 
analog  input  enters  the  DSO,  the  analog-to-digital  converter  (ADC)  selects  or  samples  the 
value  of  the  signal,  a  voltage,  at  discrete  points  in  time.  The  voltage  is  then  digitized  and 
displayed  and/or  recorded  to  be  used  later.  This  acquisition  process  generates  a  quantized 
signal,  whieh  is  a  mapping  from  a  continuous  domain  space  to  a  diserete  range  space.  By 
acquiring  many  observations  of  the  same  transient,  one  can  improve  the  S/N. 

Several  properties  of  the  DSO  affect  the  quality  of  the  data  it  can  produce.  The  ability 
of  a  DSO  to  accurately  measure  short  duration  signals  is  determined  by  the  maYtmnm 
sampling  rate  and  the  maximum  analog  bandwidth.  The  maximum  sampling  rate  indicates 
the  number  of  samples  per  second  the  ADC  can  acquire.  The  larger  the  sampling  rate,  the 
greater  the  time  resolution  the  DSO  can  display.  The  rise  time  of  the  DSO  is  a  measure  of 
the  DSOs  ability  to  respond  to  rapid  changes  in  signal.  It  is  defined  by  the  following 
formula: 
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(2.12) 


18 


where  is  the  rise  time  in  nanoseconds  and  is  the  analog  bandwidth  in  megahertz, 
which  specifies  the  frequency  range  the  DSO  can  accurately  measure.^*  Rise  times  on  good 
DSOs  with  sampling  rates  of  2  GS/s  can  be  as  short  as  800  ps. 

2.2.2  Phase-modulated  Techniques 

Phase-modulated  spectroscopy  is  a  method  of  indirectly  measuring  fluorescence 
lifetime  by  rapidly  modulating  the  excitation  light  while  simultaneously  modulating  the 
high  voltage  to  the  detector  at  almost  the  same  rate.^^  The  phase  of  the  emission  modulator 
is  shifted  relative  to  the  excitation  to  extract  the  phase  dependence  of  the  fluorescence 
intensity.  With  commercially  available  equipment  capable  of  modulation  frequencies  of 
250  MHz,  one  can  measure  lifetimes  as  short  as  1  ps. 

The  theory  and  practice  of  phase-modulated  spectroscopy  has  been  discussed  in  the 
literature.^^’^^  The  time-dependent  fluorescence  signal,  D(t),  generated  in  a 
phase-modulated  instrument  with  angular  modulation  frequency,  co  (equal  to  2%f,  where/ 
is  the  modulation  frequency  in  hertz),  and  phase-shift  angle,  (j) ,  is  given  by 

D(t)  =  A'[l-^'m^msin{cot-(j>)),  (2.13) 

where  A'  is  the  wavelength-dependent,  steady-state  (DC)  component  of  the  fluorescence 
emission,  is  the  modulation  depth  of  the  excitation  (i.e.,  the  ratio  of  the  AC  amplitude 
to  the  DC  intensity),  and  m  is  the  demodulation  factor  (given  by  the  modulation  depth  of 
the  fluorescence  divided  by  niex). 

The  fluorescence  lifetime  may  be  calculated  fi'om  the  phase  shift,  ,  or  the 
demodulation,  ,  using 
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(2-14) 
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or 
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respectively. 

The  phase-modulated  method  provides  a  fast,  precise,  and  accurate  way  to  measure 
fluorescence  lifetimes  in  the  laboratory.  Unfortunately,  when  considering  UV-visible 
fluorescence  measurements,  the  phase-modulated  method  has  lower  S/N  than  the  pulsed 
methods.  While  the  state-of-the-art  is  advancing  rapidly,  in  its  present  stage  of 
development,  phase-modulated  instrumentation  is  more  expensive  and  less  amenable  to 
fiber-optic  and  field  experiments  than  the  DSO-PMT  method  used  in  this  work. 
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3.  FACTOR  ANALYSIS 


Malinowski  has  defined  factor  analysis  as  “a  multivariate  technique  for  reducing 
matrices  of  data  to  their  lowest  dimensionality  by  the  use  of  orthogonal  factor  space  and 
transformations  that  yield  predictions  and/or  recognizable  factors.”  Some  might  consider 
the  phrase  “use  of  orthogonal  factor  space”  overly  restrictive,  since  many  factor  analysis 
methods  employ  non-oithogonal  factors;  but  this  definition  provides  an  excellent  starting 
point.  Factor  analysis  methods  for  fluorescence  data  can  be  roughly  separated  into  two 
areas:  calibration  (and  prediction)  and,  what  will  be  called  here,  profile  extraction. 

The  aim  of  calibration  is  to  extract  quantitative  information  (e.g.,  concentrations)  about 
the  constituents  of  an  unknown  sample.  Calibration  is  usually  performed  as  a  two-step 
process.  First,  in  the  calibration  step,  a  concentration-response  model  is  constructed  from 
the  instrumental  response  to  standards.  In  the  following  prediction  step,  the  instrumental 
response  of  an  unknown  is  measured,  and  the  concentration  of  the  unknown  is  estimated, 
using  the  model  from  the  calibration  step. 

The  objective  of  profile  extraction  is  qualitative  information,  in  the  form  of 
recognizable  factors,  about  the  constituents  of  an  unknown  sample.  These  recognizable 
factors  would  be  in  the  form  of  pure  component  spectra,  time  decay  profiles,  etc.,  which 
Malinowski^^  would  identify  as  real  factors.  For  example,  from  an  EEM  of  a 
two-component  solution  containing  benzene  and  naphthalene,  one  would  attempt  to  extract 
the  benzene  and  naphthalene  excitation  and  emission  spectra.  Profile  extraction  can  be 
performed  with  or  without  standards.  As  Malinowski’s  deEnition  suggests,  there  is 
sometimes  overlap  between  the  processes  of  calibration  and  profile  extraction. 
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For  environmental  applications,  one  may  have  either  or  both  goals  in  mind.  The 
concern  might  be  what  toxic  materials  are  present  and  at  what  concentrations  in  a  soil  or 
groundwater  sample.  To  accomplish  the  quantitative  goal  using  fluorescence  data,  one 
must  develop  a  calibration  model  based  on  the  relationship  between  fluorescence  intensity 
and  concentration.  The  qualitative  goal  is  accomplished  by  extracting  predictions  of 
spectra  or  fluorescence  decay  curves  and  comparing  them  with  entries  in  a  database  to 
identify  the  compound(s).  For  either  goal,  a  sound  understanding  of  the  physical  and 
mathematical  nature  of  the  fluorescence  data  is  an  important  factor  in  constructing  effective 
factor  analysis  schemes. 

The  mathematical  nature  of  different  types  of  fluorescence  data  is  addressed  in  this 
chapter.  Methods  of  calibration  and  profile  extraction  based  on  the  type  of  data  collected 
are  also  presented. 

3.1  Fluorescence  Data  as  Tensors 

Tensorial  notation  was  first  introduced  to  describe  multivariate  calibration  in  the  field 
of  chemometrics  by  Sanchez  and  Kowalski.^^’^^  Since  then,  Kowalski  et  a], have 
disseminated  tensor  terminology,  borrowed  from  mathematics  and  physics,  to  categorize 
the  types  of  data  acquired  by  various  instruments  used  by  chemists.  Objections  to  the  use 
of  the  term  tensor  in  chemometrics  have  been  raised  because  there  are  some  tensor 
conventions  which  do  not  carry  over  into  chemometrics.^®  However,  there  appear  to  be 
enough  proponents  to  carry  it  into  more  general  use.  Owing  to  Kowalski’s  influence  in  this 
area  of  chemometrics,  much  of  what  follows  is  modeled  after  papers  he  has  published  on 
this  subject.^ 
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Tensor  order  in  chemometrics  indicates  the  number  of  modes  inherent  in  a  data  set. 


Tensor  order  notation  follows  the  same  pattern  as  array  order  used  in  this  thesis  (see 
Definition  A.22).  The  notation  can  be  further  extended  to  the  instrumentation  that  collects 
the  data  to  designate  the  tensor  order  of  an  instrument. 

The  order  of  an  instrument  provides  insight  into  the  flexibility  one  has  available  in  the 
data  analysis  and  the  type  of  information  to  be  gleaned  from  the  data  produced  by  such  an 
instrument.  The  order  of  an  instrument  indicates  the  maximum  tensor  order  of  data  that 
can  be  obtained  on  a  single  sample  of  analyte.  Here,  the  focus  will  be  on  fluorescence 
instruments  and  the  various  orders  of  data  one  can  collect  with  them. 

3.1.1  Zeroth-order  Instrumentation 

Measurement  of  fluorescence  intensity  of  a  single  luminescent  species  in  solution  at  a 
single  excitation  wavelength  and  a  single  emission  wavelength  yields  a  scalar  datum,/, 
which  is  a  zeroth-order  tensor.  An  instrument  which  can  only  measure  zeroth-order  data  is 
appropriately  called  a  zeroth-order  instrument.  An  example  of  a  zeroth-order  fluorescence 
instrument  would  be  one  with  a  mercury  penlamp  excitation  source  and  a  combination 
optical  filter  and  PMT  detector.  The  concentration  of  the  luminescent  species  present  is 
proportional  to/.  Therefore,/is  a  quantitative  datum,  but  it  provides  no  qualitative 
information  concerning  the  luminescent  species.  In  fact,  the  qualitative  information,  viz., 
that  only  a  single  luminescent  species  is  present,  would  have  to  be  known  before  the 
fluorescence  measurement.  Zeroth-order  data  are  thus  very  limited  in  terms  of  the 
information  they  provide  to  the  analyst. 
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3.1.2  First-order  Instrumentation 


Measurement  of  fluorescence  intensity  as  a  tunction  of  emission  wavelength.  /(A) ,  is 
a  concept  familiar  to  chemists.  A  decade  ago,  data  acquisition  with  continuous  scanning 
monochromators  and  strip  chart  recorders  was  common.  The  practice  of  most 
spectroscopists  today  is  to  measure  spectra  at  discrete  wavelengths  with  electronic 
recording  of  the  data.  Borrowing  a  term  from  signal  processing,  we  refer  to  this  as  a 
quantized  signal.  An  ordered  collection  of  these  quantized  signals.  /(A,) .  constitutes  a 
vector,  f,  directed  from  the  origin  in  n-dimensional  space  (n-space),  i.e., 

f  _  ^{^2) 

i  ’  (3.1) 

AK)_ 

where  the  subscript  denotes  the  ith  wavelength  at  which  an  intensity  was  measured.  A 

quantized  spectrum  represents  first-order  data,  and  an  instrument  capable  of  measuring  at 
best  such  a  spectrum  is  a  first-order  instrument. 

An  excitation  spectrum,  which  is  obtained  by  recording  the  fluorescence  intensity  at  a 
fixed  emission  wavelength  while  the  excitation  wavelength  is  varied,  also  represents 
first-order  data.  Fluorescence  decay  data.  f(t) .  are  also  first-order  data  and  are  measured 
by  a  first-order  instrument. 

Two  important  properties  of  a  vector  are  its  length  and  its  direction.  When  a  single 
luminescent  species  is  present,  the  length  or  norm  of  the  spectral  vector,  ||f  || ,  provides  a 
measure  of  the  amount  of  emitter  present.  A  valuable  attribute  of  the  vector  over 
zeroth-order  data  is  that  it  provides  many  estimates  of  the  intensity  in  a  single  experiment. 
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First-order  data  therefore  provide  a  form  of  built-in  signal  averaging.  The  other  main 
property  of  the  spectral  vector,  its  direction,  gives  insight  into  the  identity  of  the  emitter; 
thus,  with  first-order  data,  one  gets  the  benefit  of  qualitative  information. 

First-order  fluorescence  data  also  include  the  emission  intensity  as  a  function  of  (1) 
excitation  wavelength,  (2)  fluorescence  decay  time,  or  (3)  polarization  of  the  emitted 
photons  relative  to  the  excitation  photons. 

3.1.3  Second-order  Instrumentation 

One  can  generate  second-order  data  by  combining  first-order  measurement  modes.  For 
example,  as  mentioned  in  the  introduction,  an  EEM  is  obtained  by  scanning  the  emission 
wavelength  at  several  different  excitation  wavelengths.  Similarly,  one  could  vary 
wavelength  and  time  to  form  a  WTM.  Such  data  can  provide  qualitative  and  quantitative 
information  like  first-order  data,  but  with  the  additional  benefits  described  in  Section  3.2.3. 
Second-order  fluorescence  instruments  capable  of  measuring  EEMs  have  been  in  use  for 
some  time.  Commercial  versions  of  these  instruments  are  now  common.  Instruments  for 
measuring  WTMs  on  a  nanosecond  or  shorter  time  scale  are  a  more  recent  innovation  and 
are  not  commercially  available  at  this  time. 

3.1.4  Third-order  Instrumentation 

A  three-mode  array  (3-array)  of  data  called  a  time-resolved  excitation-emission  matrix 
(TREEM)  combines  the  modes  of  excitation  wavelength,  emission  wavelength,  and  time 
domain.  Instruments  which  generate  these  data  are  usually  “souped-up”  versions  of  the 
second-order  instruments  which  produce  WTMs  and  WFMs,  whose  third-order  descendant 
is  the  excitation-emission  frequency  array  (EEFA).  If  time  domain  data  are  available,  the 
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modification  to  scan  both  excitation  and  emission  wavelengths  is  rather  trivial.  The 
ROST™  is  a  third-order  instrument  since  it  can  collect  TREEMs. 

The  advantages  of  this  type  of  data  will  be  discussed  in  Section  3.3.2. 

3.1.5  Fourth-order  Instrumentation 

Optical  spectroscopic  instruments  designed  to  obtain  fourth-order  data  are  possible,  but 
are  of  less  obvious  usefulness.  A  time-resolved  excitation-emission-anisotropy  instrument 
may  be  of  interest  to  biochemists^  interested  in  the  properties  of  fluorophores  in  living 
cells  and  tissue. 

3.2  Tensorial  Calibration 

One  seeks  by  calibration  to  establish  a  mathematical  and  statistical  model  of  the 
relationship  between  experimental  data  and  analyte  concentration  to  estimate 
concentrations  in  an  unknown  sample.  The  usual  procedure  is  to  create  a  set  of  standards, 
with  known  concentrations,  and  measure  the  spectral  response.  The  model  is  then  created; 
and,  in  a  prediction  step,  the  concentration  of  the  unknown  is  estimated.  A  linear 
relationship  between  the  signal  from  the  instrument  and  the  desired  parameter  is  desirable 
since  it  simplifies  data  analysis.  The  ability  to  determine  concentrations  of  multiple 
analytes  in  mixtures  and  performance  of  the  calibration  model  in  the  presence  of 
interferences  is  influenced  by  the  tensor  order  of  the  data. 

An  nth-order  instrument  generates  nth-order  data  by  definition.  Multiple  samples  of 
calibration  data  (standards)  can  be  assembled  into  an  (n  +  l)-order  tensor.  A  calibration  set 
of  individually  zeroth-order  data,  i.e.,  scalars,  is  a  first-order  tensor,  i.e.,  vector;  a  set  of 
first-order  calibration  data  set  represents  a  second-order  tensor;  etc.  However,  note  that  the 
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order  of  calibration  follows  from  the  order  of  the  instrument  and  not  from  the  order  of  the 


calibration  data  set. 

3.2.1  Zeroth-order  Calibration 

Zeroth-order  calibration  is  a  very  common  feature  in  data  analysis.  The  methodology  is 
simple,  the  analysis  is  represented  very  well  graphically,  and  it  can  easily  be  performed 
with  a  hand  calculator.  The  statistics  are  also  simple,  well-defined,  and  easy  to  calculate. 
Unfortunately,  there  is  also  very  little  flexibility  with  zeroth-order  calibration.  One  cannot 
use  zeroth-order  data  to  simultaneously  analyze  multiple  analytes  in  a  mixture.  Any 
interference  beyond  a  constant  offset  in  the  data  distorts  the  model.^’^  Interferences  present 
in  the  unknown  sample,  but  not  represented  in  the  calibration  set,  cannot  be  detected  in  the 
prediction  and  will  bias  the  parameter  (concentration)  estimate.  A  measurement  which  is 
highly  selective  for  the  target  analyte  is  the  most  favorable  case  for  zeroth-order  calibration. 

3.2.2  First-order  Calibration 

First-order  calibration  is  both  more  powerful  and  more  complex  than  zeroth-order 
calibration.  Fortunately,  there  are  several  analysis  algorithms  whose  mathematical  and 
statistical  properties  have  been  extensively  investigated.  Among  these  are  classical 
multivariate  least  squares  (MLS),  inverse  least  squares  (ILS),  principal  component 
regression  (PCR),  and  partial  least  squares  (PLS),  with  PLS,  PCR,  and  MLS  being  the  most 
popular.  First-order  data  allows  for  analysis  of  multicomponent  mixtures  with  a  maximum 
number  of  analytes  equal  to  the  number  of  elements  in  the  spectral  vector.  Interferences 
can  also  be  accommodated  in  the  calibration  model  when  using  ILS,  PLS.,  and  PCR. 

ILS,  PCR,  and  PLS  are  actually  all  inverse  least  squares  techniques,^'^’'**’'*^  in  that  they 
model  concentration  as  a  function  of  instrument  response,  opposite  to  the  way  analytical 
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expressions  are  usually  presented.  For  example,  Beer’s  Law  is  often  written 

A  =  elc,  (3.2) 

where  A  is  the  absorbance,  e  is  the  molar  absorptivity  or  extinction  coefficient,  I  is  the  path 
length  of  the  solution,  and  c  is  the  concentration.  The  inverse  relationship  may  be  written 
as 

c  =  M,  (3.3) 

where  b  is  a  model  parameter  which  defines  the  inverse  relationship.  The  advantage  of  this 
methodology  is  that  the  analysis  is  invariant  with  respect  to  the  number  of  analytes  in  the 
model.  The  concentration  of  any  spectrally  active  species  may  be  determined 
independently  of  the  others  provided  their  spectral  responses  are  modeled.  In  fact,  in  the 
method  referred  to  as  PLSl,^^  each  analyte  in  a  multicomponent  mixture  is  modeled 
separately,  and  predictions  on  the  unknown  mixture  are  conducted  independently  in  turn. 
The  analysis  for  each  component  treats  the  mixture  as  if  it  contains  only  one  analyte  and 
numerous  interferences.  The  limitation  in  this  case  is  that  each  interference  present  in  the 
unknown  must  be  included  in  the  model  or  the  prediction  will  be  biased.  However,  the 
presence  of  interferences  in  the  unknown  but  not  in  the  model  can  be  detected  by 
comparing  the  unknown  spectrum  to  the  standard  spectra. 

Data  pretreatment  is  another  area  given  much  attention  in  the  literature.  Methods  such 
as  mean  centering  and  variance  scaling  of  the  data  are  the  two  most  often  discussed.  Mean 
centering  is  accomplished  by  subtracting  the  average  value  of  the  observations  of  a  variable 
from  each  observation  for  that  variable,  i.e., 

fy  ( centered)  =  f..  -  f.  (3.4) 
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where  f.j  is  the  fluorescence  intensity  of  the  ith  sample  at  the  ^h  wavelength  and  fj  is  the 


mean  of  the  fluorescence  intensities  at  the  ^  wavelength,  i.e., 


N 
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(3.5) 


Mean  centering  accommodates  data  with  a  non-zero  baseline.  Scaling  to  unit  variance  is 
accomplished  by 


f..(scaled)  =  — 
J  s » 

J 


(3.6) 


where  Sj  is  the  standard  deviation  of  the  jth  wavelength,  i.e.. 


Performing  both  mean  centering  and  scaling  to  unit  variance  is  referred  to  as 
standardization  in  statistical  texts  and  autoscaling  to  unit  variance  in  chemometrics 
literature.  Standardizing  data  is  usually  performed  on  data  which  are  measured  on  different 
scales  with  large  range  differences  or  if  the  units  are  not  compatible. 

With  the  appropriate  data  pretreatment  and  the  optimum  algorithm,  first-order 
calibration  can  be  very  powerful.  But,  as  one  might  expect,  even  higher  order  data  increase 
versatility. 
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3.2.3  Second-order  Calibration 


The  simplest  second-order  data  follow  a  bilinear  model,  i.e.,  the  second-order  data  array 
for  a  pure  component  can  be  decomposed  into  a  dyad.  If  the  data  conform  to  the  bilinear 
model,  calibration  yields  a  third-order  tensor  in  the  calibration  model  which  is  trilinear, 
Trilinear  data  are  highly  desirable  since  they  allow  for  the  following:  (1)  calibration  with  a 
single  standard,  (2)  calibration  and  prediction  in  a  single  algorithm,  (3)  calibration  and 
accurate  prediction  in  the  presence  of  interferences,  (4)  decomposition  of  the  array  to 
provide  factor  matrices  containing  the  pure  component  spectra  of  the  two  spectral  modes, 
and  (5)  the  possibility  of  a  unique  decomposition  of  the  array.  This  kind  of  second-order 
data  is  very  powerful,  especially  in  terms  of  the  qualitative  information  that  can  be  derived. 
Unfortunately,  the  mathematical  analysis  is  complex,  and  the  statistics  are  not  well 
understood.  The  complexity  is  even  greater  if  the  second-order  data  are  not  bilinear. 

In  cases  where  pure  component  second-order  data  are  not  bilinear,  the  method  of 
analysis  will  depend  more  on  the  chemical  and  physical  nature  of  the  system.  Some 
systems  are  amenable  to  analysis  with  restricted  Tucker  models  which  may,  under  certain 
constraints,  result  in  a  unique  decomposition  and  pure  component  spectral  profiles. 
Calibration  in  other  non-bilinear  systems  by  methods  such  as  non-bilinear  rank  annihilation 
(NBRA)'*^  or  residual  bilinearization  (RBL)^’'*^  may  provide  accurate  predictions,  but  the 
ability  to  obtain  qualitative  information  is  lost  in  the  former  method. 

3.3  Profile  Extraction 

When  employing  a  profile  extraction  technique,  one  attempts  to  obtain  the  underlying 
first-order  profiles  from  the  data.  These  profiles  may  be  pure  component  emission  or 
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excitation  spectra,  time  decay  profiles,  or  relative  concentrations  in  different  samples. 
Profile  extraction  provides  the  qualitative  information  one  can  use  to  identify  an  unknown 
luminescent  species  in  a  solution. 

This  type  of  analysis  requires  data  in  the  form  of  an  iV-array,  with  N  >2,  formed  by  a 
collection  of  first-order  data,  e.g.,  first-order  calibration  data;  or  it  could  be  second-order 
data  from  a  single  solution,  e.g.,  an  EEM.  This  type  of  analysis  is  usually  labeled  by  the 
order  of  the  data  undergoing  analysis  rather  than  the  order  of  the  instrument,  which  makes 
sense  because  it  conveys  the  most  specific  information  concerning  the  capabilities  of  the 
analysis.  Thus,  profile  extraction  on  an  EEM,  which  is  a  2-array,  is  called  2-mode  profile 
extraction,  and  profile  extraction  on  a  TREEM  is  called  3-mode  profile  extraction. 

Profile  extraction  is  necessary  only  if  there  is  more  than  one  fluorescing  species  in  the 
sample.  In  cases  where  there  is  only  a  single  fluorescing  species,  the  analysis  is  trivial.  As 
in  calibration,  profile  extraction  capability  increases  with  the  order  of  the  Hata 
3.3.1  Two-mode  Profile  Extraction 

Two-mode  profile  extraction  is  the  most  commonly  applied  factor  analysis  method. 
Perhaps  this  is  so  because  there  are  so  many  ways  to  construct  first-order  data  and  many  of 
the  techniques  have  been  around  for  a  long  time.  Malinowski^^  divides  the  methods  of 
factor  analysis  into  three  categories:  (1)  abstract  factor  analysis  (AFA),  (2)  target  factor 
analysis  (TEA),  and  (3)  special  methods. 

AFA  involves  transformations  of  sets  of  abstract  factor  matrices  to  obtain  real  factors. 
Abstract  factors  are  obtained  by  eigenanalysis  of  the  data  matrix.  This  process  is  often 
called  principal  component  analysis  (PCA)  or  principal  factor  analysis  (PFA).  PFA 
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involves  the  decomposition  of  the  data  matrix,  D,  into  two  orthogonal  matrices,  R  and  C, 
often  referred  to  as  the  scores  and  loadings  matrices,  respectively,  with 

D  =  RC'.  (3.8) 

Malinowski^^  describes  a  number  of  methods  of  decomposing  D,  One  of  the  most  popular 
ones  is  singular  value  decomposition  (SVD),  in  which  D  is  represented  by 

D  =  USV'.  (3.9) 

The  vectors  in  U  are  commonly  called  the  left  singular  vectors  and  are  the  eigenvectors  of 
DD' ;  those  in  V  are  the  right  singular  vectors  and  are  the  eigenvectors  of  D'D ;  and  S  is 
the  diagonal  matrix  of  singular  values,  which  are  the  square  roots  of  the  eigenvalues  of 
DD'  (and  D'D ).  Equations  (3.8)  and  (3.9)  are  related  by  setting  R  =  US  and  C  =  V . 

The  column  vectors  of  R  and  C  represent  different  amounts  of  variance  explained  by  the 
decomposition. 

The  term  variance  in  multivariate  analysis  describes  the  amount  of  information 
contained  in  a  matrix.  It  stems  from  the  common  use  of  the  covariance  matrix  in 
multivariate  analysis.  If  G  is  the  Gramian  association  matrix  formed  from  D,  the  total 
variance  of  D  is  the  trace  of  G,  i.e.,  tr(G).  The  total  variance  of  a  matrix  is  equal  to  the  sum 
of  the  eigenvalues,  X,  of  its  decomposition,  i.e., 

tr(G)  =  2^-  (3.10) 

The  magnitude  of  the  individual  eigenvalues  represent  the  fraction  of  variance  described  by 
the  respective  eigenvectors. 

If  those  factors  that  account  for  the  bulk  of  the  variance  are  retained  and  the  others  are 
eliminated,  the  model  in  (3.8)  can  be  reduced  to  abstract  representations  of  the  data 
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resulting  from  only  real  signal,  not  from  noise.  This  objective  leads  to  attempts  to 
determine  the  array  rank  (see  Definition  A.28)  of  the  matrix.  The  term  pseudorank  is  also 
used,  since  a  data  matrix  with  experimental  noise  will  almost  always  be  full  rank  in  the 
strictest  mathematical  sense.  There  are  a  variety  of  methods  to  deduce  the  number  of 
significant  factors  in  a  decomposition  of  this  type. 

The  first  two  methods  examine  the  error  represented  in  the  eigenvalues.  Malinowski^^ 
describes  a  factor  indicator  function,  that  estimates  which  eigenvalues  of  the  data  originate 
from  error  (or  from  noise).  If  the  eigenvalues,  ,  are  ordered  from  largest  to  smallest  ( X., 
being  the  largest),  the  indicator  value  for  the  nth  eigenvector  is  given  by 


IND(n)  = 


(c-n) 


r  c 

Z^- 

j=/z+l 

ic-n) 


(3.11) 


V  J 

where  c  and  r  are  the  dimensions  of  the  data  matrix  with  c  <  r .  The  rank  of  the  data  matrix 
is  the  index  n  for  which  IND(n)  is  minimized. 

Malinowski'*^  also  developed  an  F-test  based  on  error  represented  in  eigenvalues.  The 
F-value  is  computed  from 


F(l,c-n)  = 


£(r-;  +  l)(c-7  +  l) 

;=n+l  A 


(r-n  +  l)(c-n  +  l) 


y=rt+l 


(3.12) 


The  eigenvalues  are  examined  from  smallest  to  largest.  The  index  of  the  first  eigenvalue 
whose  F(1,c  -  n)  exceeds  the  tabular  value  of  F  at  the  desired  significance  level  is  the 
estimated  rank  of  the  matrix. 
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Other  methods  examine  the  eigenvectors.  The  first-lag  autocorrelation 
coefficient  C(u) ,  which  has  been  used  to  separate  signal  eigenvectors  from  noise 

eigenvectors,  can  be  calculated  in  the  following  way: 

C 

Z“;  “w 

C(u)  =  ^!S -  (3.13) 

uu 

where  c  is  the  dimension  of  the  eigenvector,  u.  Proponents^’  of  this  method  assert  that  the 
eigenvectors  corresponding  to  real  signal  will  be  smoother  than  those  of  the  eigenvectors 
associated  with  noise.  Thus,  the  signal  eigenvectors  possess  larger  values  of  C(u) .  An 
autocorrelation  coefficient  value  of  0.5  has  been  suggested  as  a  cutoff  between  real  and 
error  eigenvectors. 

After  the  number  of  factors  has  been  determined,  the  selected  factors  are  transformed 
with  a  transformation  matrix,  T.  Transformation  is  straightforward. 

X  =  RT  (3.14) 

and 


Y'  =  T'’C'  (3.15) 

are  the  transformations  of  R  into  X  and  C  into  Y.  T  is  a  square  matrix  whose  dimensions 
equal  the  number  of  selected  factors.  Its  elements  are  usually  determined  by  iterative 
transformations  until  the  vectors  X  and  Y  meet  some  desired  criterion.  The  two  basic  types 
of  transformations  are  orthogonal  and  oblique.  Orthogonal  transformations  preserve  the 
angular  dependence  between  the  original  vectors.  Oblique  transformations  allow  variation 
of  the  angles  between  the  vectors. 
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Orthogonal  transformations  are  primarily  used  to  find  clusters  or  patterns  in  the  data  at 
hand.  An  example  would  be  trying  to  establish  a  relationship  between  chemical 
composition  and  physical  properties  of  lake  sediment.  However,  an  orthogonal 
transformation  would  not  allow  one  to  extract  physically  meaningful  spectra  from  an  EEM. 

Oblique  transformations  are  another  way  to  find  clusters  in  data  in  order  to  identify 
correlated  behavior.  Based  on  minimization  criteria  alone,  oblique  rotations  seldom 
reproduce  identifiable  spectra.  TFA^^  is  a  method  which  employs  oblique  transformations, 
although  the  transformations  matrix  is  determined  in  an  different  way.  Also,  a  constraint 
included  in  the  transformation  algorithm  can  produce  good  results. 

TEA  attempts  to  extract  meaningful  factors  one  at  a  time  by  comparing  transformations 
of  the  abstract  factors  to  target  factors,  i.e.,  real  factors  one  suspects  are  part  of  the  data.  In 
TEA,  one  generates  a  transformation  vector,  t ,  by  fitting  the  target  vector,  x ,  to  the 
abstract  factor  matrix  in  a  least  squares  sense,  i.e., 

t  =  R^x.  (3.16) 

The  predicted  vector,  x ,  is  then  obtained  by 

x  =  Rt.  (3.17) 

Equations  (3.16)  and  (3.17)  can  be  combined  as 

x  =  RR‘^x  (3.18) 

to  reveal  that  x  is  a  least  squares  projection  of  x  into  the  space  spanned  by  the  data  as 
represented  by  the  abstract  factors.  The  predicted  vector  is  subjected  to  a  statistical  test  (F- 
test)  to  determine  whether  it  is  a  valid  real  factor.  TEA  requires  that  the  experimenter  have 
some  idea  of  what  may  be  in  the  sample  and  to  have  target  vectors  on  hand.  Of  course,  if 
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one  has  a  database,  this  should  not  be  a  problem,  but  for  an  unknown  which  may  have 
spectra  not  in  the  database,  this  would  pose  a  problem. 

Special  methods  of  factor  analysis  encompass  a  wide  variety  of  techniques.  The  most 
common  are  classified  as  evolutionary  factor  analysis  (EFA)  methods,  which  include 
modeling  and  self-modeling  algorithms.  EFA  methods  often  employ  constraints  based  on 
known  physical  properties  of  the  measured  parameters.  A  common  constraint  for 
wavelength-mode  spectroscopic  data  is  nonnegativity Exponential  decay  is  the  usual 
constraint,  or  model,  for  time  mode  data  in  spectroscopy  and  kinetics.^^'^^  Unimodality  and 
nonnegativity  have  been  used  as  constraints  in  chromatography.  Modeling  and 
self-modeling  methods  of  EFA  lend  themselves  readily  to  analysis  of  EEMs  and  WTMs. 

Modeling  methods  can  be  used  very  effectively  with  fluorescence  lifetime  data.  One 
simply  uses  the  model  in  (2.1 1)  to  fit  the  time  domain  data.  This  technique  was  used  by 
Knorr  and  Harris^^  to  extract  lifetime  and  emission  spectra  fi:om  WTMs  of  two-component 
solutions.  Fluorescence  decay  following  pulsed  excitation  was  measured  at  a  number  of 
wavelengths.  The  WTM,  D,  is  represented  as  a  bilinear  array, 

D  =  WT',  (3.19) 

where  W  is  the  matrix  of  wavelength  mode  factors  and  T  is  the  matrix  of  time  mode 
factors.  The  convoluted  decays  were  fit  in  a  Simplex  search  algorithm  using  a  single 
lifetime  parameter  for  each  fluorescing  species.  Following  each  iteration  of  the  simplex 
algorithm,  Knorr  and  Harris  obtained  a  least  squares  fit  of  the  time  mode  vector  to  estimate 
the  wavelength  factor  matrix,  i.e., 

W  =  DT(T'T)'’ .  (3.20) 
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Many  practitioners^'^’^^  now  call  this  technique  global  analysis.  It  is  finding  extensive  use 
in  biochemistry,  especially  kinetics.  Often,  users  will  search  parameter  space  for 
parameters  in  the  time  mode  and  the  wavelength  mode  while  employing  the  Marquardt^^ 
search  algorithm,  although  a  faster  scheme  would  employ  the  Marquardt  method  and  the 
least  squares  method  as  Knorr  and  Harris  did.  This  allows  a  faster  search  of  parameter 
space  for  the  nonlinear  part  (lifetime)  and  a  fast  solution  to  the  linear  part  (wavelength).  Of 
course,  if  both  modes  are  built  on  nonlinear  models,  such  a  scheme  would  not  be  possible. 
One  could  also  use  nonnegativity  constraints  on  the  values  of  the  wavelength  mode  to  force 
a  realistic  solution. 

Self-modeling  methods  impose  constraints  on  the  analysis  without  using  any  strict 
model  for  the  behavior  of  the  data.  For  example,  constraining  an  oblique  transformation  of 
abstract  factors  to  produce  only  positive  values  in  the  wavelength  domain  has  been  used  to 
decompose  fluorescence  of  several  mixtures  into  recognizable  spectra.^^’^® 

Self-modeling  curve  resolution  is  a  popular  form  of  this  type  of  analysis.  Setting 
constraints  such  as  nonnegativity  for  spectra  and  using  an  alternating  least  squares 
algorithm  have  been  used  with  success.'*®'^^  As  usual,  the  data  matrix,  D,  is  represented  as 
the  product  of  two  factor  matrices,  A  and  B,  i.e., 

D  =  AB' .  (3.21) 

With  some  starting  value  for  A,  B  is  estimated  by  least  squares,  i.e., 

B'  =  A^D,  (3.22) 

and  any  elements  of  B  that  are  negative  are  set  to  zero.  Then,  A  is  estimated  from  B  using 


and  any  elements  of  A  that  are  negative  are  set  to  zero.  The  process  is  continued  until  a 
termination  criterion,  such  as  minimum  residual  (see  Appendix  B),  is  reached.  One  could 
also  employ  a  nonnegative  least  squares  (NNLS)  algorithm. 

EFA  procedures  are  simple  to  understand,  and  they  do  not  require  prior  knowledge 
about  the  components  in  the  system.  By  following  a  simple  model,  e.g.,  exponential  decay 
or  imposing  constraints  and  iterating,  one  can  generate  the  real  factors  in  many  cases. 
However,  there  is  a  danger  that  one  could  get  a  result  that  is  not  the  correct  outcome  as  a 
consequence  of  what  is  called  the  rotation  problem. 

The  rotation  (or  transformation)  problem  is  a  direct  result  of  (3.14)  and  (3.15).  When 
operating  on  a  2-array,  one  can  construct  any  transformation  matrix  and  use  it  in  (3.14)  and 
(3.15).  The  result  is  that  the  decomposition  of  D  is  not  unique.  Even  employing 
constraints  such  as  nonnegativity  does  not  ensure  uniqueness.*®  However,  decomposition 
of  trilinear  3-arrays  does  provide  a  unique  solution. 

3.3.2  Three-mode  Profile  Extraction 

Third-order  data  that  follow  the  trilinear  model  can  be  decomposed  uniquely.  Some  of 
the  algorithms  for  performing  the  decomposition  employ  eigenanalysis  based  procedures 
(EBP),  while  others  use  the  parallel  factors  (PARAFAC)  analysis.  PARAFAC  is  also 
called  alternating  least  squares  (ALS)  or,  as  we  shall  henceforth  refer  to  it,  three-mode 
alternating  least  squares  (3M-ALS). 

The  formal  treatment  of  EBPs  can  be  found  in  Appendix  C.  They  are  used  to  obtain  a 
direct  solution  to  the  trilinear  decomposition.  There  are  variants  of  the  EBP,  but  these  deal 
mostly  with  the  method  of  generating  the  data  slices  to  be  used  in  subsequent  eigenanalysis. 
These  methods  are  fast  and  produce  good  results  on  synthetic  data  and  on  real  data  with 
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small  numbers  of  components  and  high  S/N.  This  methodology  has  been  applied  to  several 
types  of  spectroscopic  data.*°’^^’^^’^^'^^  However,  studies  have  shown  that  the  presence  of 
noise  can  be  problematic  and  that  EBPs  followed  by  3M-ALS  is  necessary.*®’®^ 

3M-ALS  is  an  iterative  approach  to  trilinear  analysis.  Its  methodology  is  described  in 
Appendix  B.  3M-ALS  usually  produces  good  results  even  in  the  presence  of  noise. 
Unfortunately,  the  iterative  approach  can  be  very  time  consuming,  requiring  tens  of 
thousands  of  iterations.  In  addition,  there  are  some  mathematical  curiosities  associated 
with  3M-ALS. 

Kruskal  et  al.^  reported  a  phenomenon  they  referred  to  as  two-factor  degeneracy.  A 
two-factor  degeneracy  exists  when  two  factors  in  a  trUinear  decomposition  (1)  are  highly 
correlated  in  all  three  modes,  (2)  have  very  large  and  almost  equal  magnitudes,  (3)  have  the 
overall  effect  of  cancellation  of  one  factor’s  contribution  by  the  other. 

Another  phenomenon  that  has  been  associated  with  3M-ALS  is  multiple  local  optima 
(MLO).*°  MLO  appear  when  two  starting  sets  of  parameters  yield  different  solutions  to  the 
3M-ALS  sequence.  Mitchell^®  has  suggested  that  MLO  are  a  function  of  termination 
criteria  in  many  cases. 

One  way  to  expedite  the  convergence  of  the  3M-ALS  algorithm  with  real  spectral  data 
would  be  to  combine  it  with  some  of  the  techniques  used  in  self-modeling  curve  resolution. 
By  imposing  a  nonnegativity  constraint,  the  algorithm  may  be  forced  to  converge  sooner. 
Three-mode  nonnegative  alternating  least  squares  (3M-NNALS)  may  also  help  to  resolve 
factors  difficult  to  separate  with  3M-ALS.  Such  an  algorithm  can  easily  be  created  by 
setting  all  negative  elements  in  each  mode  to  zero  after  obtaining  a  solution  for  fliat  mode.^^ 
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The  danger  is  that  the  procedure  may  not  converge  to  the  unique  solution  if  the  data  are 
noisy.  Another  problem  with  3M-ALS  is  that  there  are  few  promising  rank  estimation 
procedures  available. 

It  is  ironic  that  while  there  are  an  abundance  of  rank  estimation  procedures  for 
second-order  data  (specifically,  2-arrays),  there  is  no  mathematical  formalism  for 
determining  the  rank  of  a  3-array.^^  Thus,  rank  estimation  is  veiy  difficult  with  trilinear 
data. 

Often,  TLD  is  performed  on  a  3-array  for  a  number  of  trial  ranks.  The  rank  of  trilinear 
data  is  determined  by  examining  the  resultant  factors  for  each  trial  rank.  When  the  true 
rank  of  the  data  is  exceeded,  resolutions  should  produce  factors  that  model  noise  and  will 
appear  as  noisy  factors  in  all  three  modes. 
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4.  PATTERN  RECOGNITION 


From  the  time  people  are  very  young,  they  are  very  adept  at  pattern  recognition.  A  baby 
easily  recognizes  his  mother’s  voice  and  face.  The  phenomenal  capabilities  of  the  brain 
allow  human  beings  to  detect  extremely  subtle  differences  in  sounds,  images,  odors,  and 
textures,  to  name  a  few. 

For  data  analysis,  it  would  certainly  be  beneficial  to  train  a  computer  to  recognize 
patterns  and  to  place  samples  into  different  groups  or  classes.  Unfortunately,  this  is  a 
nontrivial  task.  What  humans  perform  innately  must  be  mathematically  modeled  for  the 
computer  to  perform.  Important  characteristics  about  the  items  to  be  grouped  must  be 
identified  and  given  a  numerical  value.  These  are  called  descriptors  of  the  data.  Also,  the 
method  one  uses  to  group  the  data  will  depend  on  what  information  one  has  beforehand. 

Pattern  recognition  methods  fall  into  two  categories:  supervised  classification  and 
unsupervised  classification.^  Supervised  classification,  or  simply  classification,  seeks  to 
place  items  into  one  of  two  or  more  known  groups.  Unsupervised  classification,  or 
clustering,  seeks  to  place  items  into  groups  based  on  similarities  without  association  to 
previously  established  groups.  To  avoid  confusion,  henceforth,  we  will  refer  to  supervised 
classification  as  classification  and  unsupervised  classification  as  clustering. 

4.1  Classification 

Classification  algorithms  are  developed  based  on  specific  rules  for  selecting 
membership  in  a  class.  While  there  are  various  rules  for  determining  class  membership,  the 
ultimate  goal  is  to  place  the  unknown  into  the  correct  class.  The  quality  of  a  classification 
rule  is  determined  by  the  probability  of  misclassification.  The  occurrence  of  noise  makes 
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it  unlikely  that  100  percent  classification  accuracy  can  be  obtained  for  experimental  data. 
When  the  probability  of  misclassification  is  small,  it  is  called  an  optimal  classification  rule, 
4.1.1  Optimal  Classification  -  Bayes’  Rule 

Several  authors  consider  Bayes’  rule  to  be  the  optimal  classification  rule.^^’^®  Bayes’ 
rule  holds  that  one  should  place  an  item,  whose  data  are  in  x,  in  the  group  with  the  highest 
posterior  probabiUty,  P(G,  |x)  (read  as  the  probability  of  belonging  to  group  G,-  given  the 


data  in  x.)  Thus,  one  wants  to  assign  the  item  to  the  group  i  such  that 

P(G,.|x)>p(g^.|x)  for  ally  (4.1) 

P(G,|x)  is  not  easily  obtained  in  practice,  though  it  may  be  realized  through  Bayes’ 
theorem,  viz.. 


P(«|C,)P(g,) 

2p(*|g,)p(g,)’ 


(4.2) 


where  P(G,. )  is  the  prior  probability  of  occurrence,  i.e.,  the  probability  of  an  item  being  a 
member  of  a  group  in  a  population  (without  regard  to  the  descriptors  used  to  classify  it), 
P(*|G,)  is  the  probabihty  of  acquiring  the  descriptors  in  x  given  the  group  G,- .  Thus,  the 


item  with  the  largest  value  for  p(x|G.  )P(Gi )  will  have  the  largest  posterior  probability,  and 
the  item  would  be  assigned  to  the  group  which  satisfies 

p(x|G,)p(G,)>p(x|G^)p(G,)  forall;#!.  (4.3) 

Assume  the  descriptors  in  the  groups  are  normally  distributed  with  mean  |x,  and 
covariance  Z,. .  Then, 
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(4.4) 


)  2:7 ‘ (x  -  H.- )] . 

where  p  is  the  number  of  descriptors  in  x  and  |Z,  |  is  the  determinant  of  Z,. .  It  can  be 
shown  that  one  can  generate  the  quadratic  discriminant  score,  dP ,  by  combining  the  left 
side  of  (4.3)  and  (4.4)  and  eliminating  constants  that  are  the  same  for  all  populations,®^  i.e., 

=-iln|Z,|-|^-|(x.n,)'2:r'(’:-t‘,)]+l“P(G,)-  (4-5) 

From  (4.3),  then,  the  item  with  descriptors,  x,  would  be  allocated  to  the  group  with  the 
largest  quadratic  discriminant  score. 

Further  assume  that  the  groups  share  a  common  covariance,  Z.  It  also  can  be  shown 
that  one  can  generate  a  linear  discriminant  score  dj  analogous  to  (4.5).®^  However,  a  more 
convenient  form  for  a  linear  discriminant  is  the  Mahalanobis  distance  squared,  , 

Sf  =(x-Iii)  (4-6) 

Allocation  of  the  item  with  x  to  group  G,  is  based  on 

--5f  +  lnP(G,)>--55  +  lnP(Gj  foralljVi.  (4.7) 

Note  that  if  the  prior  probabilities  are  equal,  the  inequality  in  (4.7)  reduces  to  allocation 
based  on  the  minimum  8f .  , 

In  practice,  the  population  covariance  (for  the  case  of  common  covariances),  Z,  and 
means,  |X, ,  are  not  usually  known.  Instead,  they  are  estimated  from  calibration  sets  of  data 
giving  the  pooled  covariance  matrix,  Sp^,^,  and  the  sample  means,  x,. ,  respectively. 

These  are  calculated  by 
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(4.8) 


c  =-i^l _ 

•spooled  £  ’ 

i=I 

where  L  is  the  number  of  groups,  rij  is  the  number  of  items  in  the  ith  group.  S,.  is  the 
sample  covariance  of  the  ith  group,  given  by 

-x,  )(x^.  -X,.) 

- 7 - - ,  (4.9) 

and  X,.  is  the  sample  mean  vector  of  the  descriptors  of  the  ^h  item  in  the  ith  group,  given 
by 

X,  =— |^x«.  (4.10) 

rii  " 

Optimal  classification  rules  may  include  additional  considerations.®^  One  of  these  is 
the  “cost”  of  misclassification,  i.e.,  a  number  representing  the  potential  value  of  a 
misclassification;  e.g.,  an  engineer  classifying  a  failing  bridge  as  safe  would  be  more  costly 
than  the  reverse. 

4.1.2  Fisher’s  Method  of  Linear  Discrimination 

One  objective  of  pattern  recognition  schemes  is  to  find  those  aspects  of  the  data  which 
best  describe  the  group.  Feature  selection  (also  referred  to  as  feature  extraction)  endeavors 
to  map  the  original  descriptors  into  a  lower  dimensional  space  while  providing  greater 
group  separation.  Canonical  discrimination  or  Fisher’s  linear  discriminant  function,®^ 
which  is  a  feature  selection  technique  for  discriminating  among  several  groups,  is  described 
in  capsule  form  here. 
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Assume  a  common  covariance  matrix  for  the  different  groups  and  sample  estimates 


rather  than  population  statistics.  Define  the  overall  mean  vector  x  as 


L  n, 


X  =  *=i 


I", 


(4.11) 


i=l 


Define  the  sample  between  groups  matrix  as 


(4.12) 


This  is  a  measure  of  the  spread  of  the  means  of  the  groups.  Next,  define  the  sample  within 
groups  matrix  W  as 

WrXSK  -X/)(Xib  -x,)  .  (4.13) 

1=1  t=i 

W  is  a  measure  of  the  spread  of  the  values  within  a  group  and  is  related  to  Spooled  ^ 
constant. 

The  linear  discriminant  functions  are  determined  by  finding  the  eigenvectors  e  of 
This  optimizes  tr(W"‘Bj),ameasure  of  class  separability.^®  This  is 
accomplished  by  solving  the  eigenvalue  equation, 

W-’Boe  =  Xe.  (4.14) 

Since  B^  has  rank,  L  - 1 ,  there  will  be  only  r  =  min(p,  L— l)  nonzero  eigenvalues,  X . 

Thus,  only  the  r  corresponding  eigenvectors  need  to  be  utilized  to  generate  features  or 
linear  discriminants  without  losing  class  separation.  Following  eigenanalysis,  the 
eigenvectors  are  each  scaled  by 
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ec  =  l, 


(4.15) 


such  that 

^  Spooled  1  ~  1  > 

providing  the  coefficients  1  for  a  linear  discriminant  function. 

The  inner  product  of  1  and  x  produces  the  sample  linear  discriminant,  y.  Using  the 
coefficients  corresponding  to  the  largest  eigenvalue,  one  generates  the  sample  first 
discriminant,  }’,  =  I'x .  With  the  p  by  r  matrix  of  the  scaled  significant  eigenvectors  L,  one 
maps  X  into  the  new  discriminant  space,  forming  the  sample  discriminant  vector  y  =  L'x . 

If  r  <  3 ,  the  sample  discriminants  may  be  plotted  against  one  another  for  each  sample, 
showing  the  groups  optimally  separated  in  the  discriminant  space. 

Using  all  of  the  r  significant  eigenvectors  provides  for  a  classification  scheme  identical 
to  (4.6)  and  (4.7).  The  Mahalanobis  distance  squared,  ,  is  given  by 

=(y-y.)  (y-yi)=[L'(x-x,)]  [l'(x-xJ].  (4.i6) 

4.2  Cluster  Analysis 

Cluster  analysis  is  a  method  of  separating  items  into  groups  based  on  where  the 
descriptors  of  the  data  place  them  in  space.  As  with  all  of  the  techniques  discussed  so  far, 
there  are  a  multitude  of  clustering  techniques.  One  fast  method  which  is  relatively  easy  to 
use  and  understand  and  provides  an  easily  interpreted  graphic  is  agglomerative  hierarchical 
cluster  analysis  (HCA).^^’’° 

HCA  operates  by  finding  similarities  (or  dissimilarities  depending  on  the  measure  one 
chooses  to  utilize)  between  items  one  at  a  time.  It  is  an  iterative  process,  taking  N-1  steps 
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for  N  items.  Groups  are  formed  as  the  items  become  linked  together  based  on  similarity. 
Each  time  a  group  is  formed,  it  is  treated  as  a  single  item.  The  measure  of  similarity  at  the 
time  each  group  is  formed  becomes  the  level  of  that  group  in  the  hierarchy  of  groups. 

Thus,  at  the  start  of  the  process,  there  are  N  groups;  and,  at  the  end,  there  is  only  one,  and 
the  levels  indicate  the  natural  groupings.  Each  time  a  new  group  is  formed,  i.e.,  a  linkage  is 
formed,  the  similarities  change  to  reflect  the  relationship  of  that  group  to  the  other  groups. 

Three  common  ways  of  forming  linkages  between  groups  are  (1)  single  linkage, 
choosing  the  two  items  within  the  different  groups  with  the  greatest  similarity  between 
them  and  assigning  that  as  the  new  measure;  (2)  complete  linkage,  choosing  the  two  items 
within  the  different  groups  with  the  least  similarity;  (3)  average  linkage,  creating  an 
average  over  all  the  items  between  the  groups.  The  choice  of  method  will  depend  on  the 
data;  however,  examination  of  all  three  methods  is  recommended.  Different  measures  of 
similarity  should  also  be  explored. 

The  best  measure  of  similarity  will  depend  on  the  nature  of  the  data  under  examination. 
A  contmion  measure  of  dissimilarity,  perhaps  useful  as  a  first  guide  for  all  data,  is  distance. 

It  is  natural  to  associate  objects  based  on  how  close  they  are  in  space.  After  all,  that  is  the 
nature  of  groups.  The  distance  between  any  two  objects  in  coordinate  space  defined  by 
their  descriptor  vectors,  d(x,y) ,  can  be  determined  using  a  variety  of  metrics.  A  versatile 
one  is  the  Minkowski  metric,®^’^°  i.e., 

d{x,y)  = 

where/?  is  the  number  of  descriptors  in  x  and  y  and  m  defines  the  specific  metric.  For 


X; 


i=l 


/m 


(4.17) 
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example,  if  m  =  1 ,  <i(x,y)  is  the  “city-block”  or  “Manhattan”  distance,  so  called  because  it 

is  analogous  to  measuring  distance  on  a  city  grid  while  driving  or  walking.  Using  m  =  2 
results  in  the  familiar  Euclidean  or  “straight-line”  distance.  Relating  each  object  to  every 
other  object  is  a  simple  task  from  the  matrix  of  distances. 

One  measure  of  similarity  which  has  been  used  for  spectra  is  the  angle  cosine  ’  (see 
Definition  A.20).  The  angle  cosine  can  take  on  values  from  -1  to  1,  although  with  real 
nonnegative  spectra,  its  values  will  only  fall  from  1  to  0.  A  high  value  of  the  angle  cosine 
indicates  more  similarity:  perfect  agreement  is  indicated  by  the  value  1.  One  notes  a 
similarity  to  the  familiar  correlation  coefficient,  and  the  two  are  related.  To  make  the 
similarity  measure  a  minimum  like  distance ,  i.e.,  small  value  equals  high  similarity,  the 
angle  cosine  is  simply  subtracted  from  one.  The  angle  cosine  has  the  added  benefit  that  the 
magnitude  of  the  individual  vectors  do  not  affect  the  correlation. 
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5.  RESULTS  AND  DISCUSSION 


The  objective  of  this  research  was  to  examine  chemometric  methods  that  can  be  applied 
to  determine  the  concentrations  of  luminescent  species  in  solutions  and  mixtures  and  the 
identity  of  composite  mixtures  from  fluorescence  data  obtained  by  the  ROST™  or  a  similar 
instrument.  Specifically,  we  investigated  the  chemometric  analysis  of  time-resolved 
multimode  fluorescence  data.  The  background  literature  in  this  area  is  almost  nonexistent, 
particularly  for  analysis  of  actual  data,  and  we  felt  that  investigation  of  a  wide  range  of 
approaches  would  yield  greater  progress  than  in-depth  studies  of  narrower  topics.  Our 
research  should  be  viewed  more  as  establishing  a  path  for  future  research  and  validation  of 
chemometric  approaches  than  as  a  detailed  comprehensive  study. 

We  shall  describe  and  discuss  chemometric  analysis  of  three  data  sets.  The  first  two 
data  sets  were  acquired  for  dilute  solutions  with  only  a  few  components.  The  first  data  set. 
Data  Set  One,  was  measured  as  a  WTM-Concentration  Array.  The  second  set.  Data  Set 
Two,  was  measured  as  a  TREEM.  Our  analysis  goals  in  each  case  were  to  extract  the 
spectral  and  temporal  profiles  for  each  luminescent  species  present  in  the  solution. 

Data  Set  Three  consists  of  WTMs  of  four  different  fuel  products  present  at  three 
different  concentrations  on  three  different  soil  matrices.  Our  analysis  goals  for  these  data 
were  to  extract  profiles  for  the  luminescent  species  in  the  fuels,  to  classify  the  fuels  using  a 
priori  information  on  the  fuel  types,  and  to  perform  cluster  analysis  to  group  the  fuels 
without  the  benefit  of  a  priori  information. 
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5.1  Profile  Extraction  of  Simple  Systems 

We  defined  a  chemically  simple  system  as  one  in  which  a  small  number  of  luminescent 
species  are  present  at  such  sufficiently  low  concentrations  that  energy  transfer  effects  are 
negligible.  Such  simple  systems  offer  the  best  prospects  for  extracting  chemically 
recognizable  factors  from  the  raw  data.  We  made  this  assumption  because  simulated  data 
utilized  in  most  tests  of  chemometric  methods  resemble  data  from  simple  systems.*®’^®’^^ 

5.1.1  Analysis  of  Three-mode  Data  With  an  EBP,  3M-ALS,  and  3M-NNALS 
We  focus  on  several  key  aspects  of  trilinear  analysis  in  this  section.  Trilinear 

three-mode  data  have  the  important  inherent  property  that  they  can  yield  rotationally  unique 
decompositions.*^  In  the  analysis  by  an  EBP  or  3M-ALS,  one  can  either  estimate  the  rank 
of  the  3-array  before  analysis  or  perform  the  analysis  for  various  assumed  array  ranks  until 
a  satisfactory  decomposition  is  obtained;  the  latter  is  the  usual  procedure. 

5.1.1.1  Data  Set  One:  A  WTM-Concentration  Array 

The  first  sample  evaluated  contained  fluorene,  pyrene,  and  naphthalene  in  aqueous 
solution.  These  compounds  were  selected  because  their  spectra  exhibit  only  modest 
overlap  along  the  wavelength  mode  and  their  hfetimes  are  well-separated.  Solutions  were 
prepared  by  successive  additions  of  aliquots  of  stock  solutions  of  the  analytes  to  3  mL  of 
water  in  a  standard  1-cm  internal  path  fused  silica  cuvette.  The  stock  solutions  were 
prepared  in  spectral  grade  methanol  at  concentrations  of  1  ppm  (±  7%),  100  ppm  (±  5%), 

and  100  ppm  (±  5%),  for  fluorene,  pyrene,  and  naphthalene,  respectively.  A  total  of  seven 
solutions  (including  the  water  blank)  were  prepared  according  to  the  design  in  Table  5.1. 
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Contour  plots  of  the  WTMs  that  were  measured  for  each  of  the  seven  solutions  are 
presented  in  Figure  5.1.  Note  that  fiber  optics  were  not  used  for  light  delivery  and 
collection  in  these  experiments.  The  excitation  source  was  the  fourth  harmonic  (266  nm) 
of  a  Nd:  YAG  laser  (Spectra  Physic  GCR-12).  The  fluorescence  emission  was  collected  at 
right  angles  to  the  excitation  and  imaged  with  a  lens  onto  the  entrance  slit  of  a  320  mm 
focal  length  monochromator  (ISA  HR-320).  Decay  profiles  were  quantized  at  2  ns 
intervals  over  a  400  ns  range;  and  emission  spectral  profiles  were  coDected  at  10  nm 
intervals  over  the  wavelength  range  280  nm  -  480  nm.  Therefore,  each  WTM  consists  of 
201x21  or  4221  data  elements. 


Table  5.1.  Experiment  matrix  for  Dal 

ta  Set  One 

SOLUTION 

Total  Fluorene  (ppb) 

Total  Naphthalene  (ppb) 

1 

0 

0 

0 

2 

0.33 

0 

0 

3 

0.33 

33 

0 

4 

0.33 

33 

33 

5 

1.0 

33 

33 

6 

1.0 

67 

33 

7 

1.0 

67 

100 

5.1.1.1.1  Rank  Estimation  of  Data  Set  One 

With  the  knowledge  that  the  solutions  were  prepared  by  standard  additions,  it  is  easy  to 
deduce  that  there  are  at  least  four  components  represented  in  the  data.  The  water  Raman 
scatter  is  the  dominant  feature  for  WTM  1  (that  of  the  blank)  in  Figure  5.1,  and  the 
additions  of  analyte  are  evident  in  WTMs  2, 3,  and  4.  The  qualitative  differences  between 
WTM  4  and  WTMs  5-7  are  not  nearly  as  easy  to  discern.  A  significant  challenge  for  the 
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WTMl 


WTM6 


WTM3 


Figure  5.1.  Contour  plots  of  seven  WTMs  of  Data  Set  One.  Contours  represent  levels  of 
equal  fluorescence  intensity.  The  weakest  intensity  is  represented  by  the  outermost  contour 
in  each  WTM;  it  has  a  value  of  0.025  (arbitrary  units).  Intensity  increases  by  0.025  for  each 
successive  contour. 
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chemometric  methods,  therefore,  is  to  accurately  predict  the  data  ranks  for  individual 
WTMs  and  for  the  entire  data  set. 

S  VD  were  tested  on  this  data  set.  Since  the  SVD  can  be  performed  for  each  individual 
WTM,  the  rank  estimation  methods  were  applied  to  each  individual  WTM.  Malinowski's 
factor  indicator  function  and  F-test  use  only  the  eigenvalues,  whereas  the  eigenvectors 
serve  as  input  if  the  rank  is  estimated  from  the  first-lag  autocorrelation  coefficient.  Note 
that  first-lag  autocorrelation  coefficients  can  be  calculated  in  both  the  wavelength-mode 
space  and  the  time-mode  space.  An  autocorrelation  value  of  0.5,  suggested  for  a  cutoff 
between  real  and  error  eigenvectors,  was  employed  here. 

The  solutions  in  Data  Set  One  were  constructed  to  contain  one,  two,  or  three  added 
chemical  components,  and  the  water  Raman  signal  represents  a  fourth  component.  The  total 
number  of  components,  which  we  refer  to  as  the  a  priori  (correct)  data  rank,  is  shown  in 
Table  5.2  along  with  the  rank  estimations  by  the  various  methods.  By  inspection  of  Table 
5.2,  it  is  clear  that  the  autocorrelation  coefficient  along  the  time  mode  is  a  totally  unreliable 
rank  predictor  at  the  recommended  0.5  cutoff  value.  The  agreement  between  a  priori  and 
estimated  rank  improves  substantially  when  the  cutoff  value  is  increased  to  a  higher  value, 
e.g.,  4-0.9;  but  at  this  stage,  there  is  no  strong  justification  for  choosing  a  different  cutoff 
value.  The  factor-indicator  method  and  the  F-test,  which  are  based  on  eigenvalues, 
accurately  predict  the  rank  for  WTMs  3  through  7,  but  overestimate  the  rank  for  WTMs  1 
and  2.  In  fact,  both  methods  actually  give  a  higher  rank  for  WTMs  1  and  2  than  for  the 
other  WTMs  even  though  the  number  of  chemical  components  is  definitely  less.  The 
wavelength-mode  autocorrelation  estimate  exhibits  reverse  behavior.  It  agrees  with  the  a 
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priori  rank  for  WTMs  1  and  2,  but  underestimates  the  rank  by  one  for  all  of  the  solutions 
that  contain  four  components.  We  speculated  that  the  fact  that  the  Raman  scatter  spectrum 
consists  of  a  single  non-zero  intensity  value  at  the  wavelength  interval  used  in  the  data 
collection  is  the  reason. 


Table  5.2.  Rank  estimates  for  Data  Set  One 
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^  Rank  is  equal  to  the  number  of  emitting  species;  parenthetical  value  includes  Raman. 
^  Malinowski's  factor  indicator. 

^  Malinowski's  F-test  for  significance  above  the  5%  level. 

^  Autocorrelation  coefficient,  values  above  +0.5. 


Unfortunately,  each  of  these  methods,  as  applied  normally,  estimates  the  rank  for  a 
2-array.  This  could  cause  a  problem  if  a  luminescent  species  is  present  in  one  solution,  but 
not  another;  then  the  rank  estimate  of  the  3-array  would  be  deficient.  For  example, 
consider  the  estimate  of  the  rank  of  a  3-array  consisting  of  two  WTMs  where  each  WTM 
contained  fluorescence  from  five  luminescent  species,  but  there  were  only  four  species 
common  to  both  WTMs.  In  such  a  case,  the  rank  of  each  WTM  could  be  correctly 
estimated  as  four,  but  the  rank  of  the  3-array  would  be  underestimated.  We  have  therefore 
tested  the  2-array  rank  estimation  methods  to  a  compressed  3-array  obtained  by  adding  all 
of  the  individual  WTMs.  The  results  are  contained  in  the  last  row  of  Table  5.2.  Note  that 
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the  method  cannot  provide  a  correct  rank  estimate  if  the  rank  of  the  3-array  is  lower  than 
the  smaller  dimension  of  the  matrix.  This  is  not  a  problem  here  since  the  smaller  dimension 
of  the  new  matrix  obtained  by  summing  the  WTMs  is  2 1 . 

This  WTM  summation  procedure  appears  to  cause  the  factor  indicator  function  and  the 
F-test  to  overestimate  the  rank  by  one.  The  wavelength-mode  autocorrelation  estimate  is 
equal  to  the  number  of  components  in  solution  plus  Raman  scatter,  but  this  may  be 
fortuitous.  Once  again,  the  time-mode  autocorrelation  is  poor  using  the  +0.5  cutoff. 
5.1.1.1.2  TLD  of  Data  Set  One 

We  next  tested  TLD  of  Data  Set  One  for  assumed  ranks  of  two  through  six.  Only 
selected  results  are  shown  here.  Note  that  TLD  treats  the  entire  data  set,  which  consists  in 
this  case  of  the  seven  WTMs  shown  in  Figure  5.1.  The  first  try  at  decomposition  used  the 
EBP  of  Sanchez  and  Kowalski,^®  and  the  rank  four  decomposition  results  are  presented  in 
Figure  5.2  along  the  emission  wavelength,  fluorescence  decay  time,  and  concentration 
modes.  Factors  Two  and  Four  have  negative  intensity  components  in  the  wavelength  and 
time  modes,  which  are  physically  impossible;  but  the  factors  can  nonetheless  be  reasonably 
assigned  as  follows:  Factor  One  -  naphthalene;  Factor  Two  -  fluorene;  Factor  Three  - 
pyrene;  Factor  Four  -  water  Raman  scatter. 

Next  we  tried  3M-ALS  on  these  same  data  using  the  factors  from  the  EBP  as  starting 
vectors.  The  convergence  criterion  for  the  3M-ALS  analysis  sequence  was  a  difference  in 
consecutive  norm  of  residuals  less  than  1  x  (this  criterion  was  used  for  all  calculations 
in  this  thesis).  The  algorithm  provided  a  converged  solution  in  104  iterations.  The  factors 
for  the  rank  four  3M-ALS  decomposition  are  displayed  in  Figure  5.3.  They  improve  upon 
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Figure  5.2.  Rank  4  TLD  of  Data  Set  One  using  EBP  of  Sanchez  and  Kowalski. 
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Figure  5.2.  Continued. 
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Figure  5.3.  Rank  4  3M-ALS  determined  factors  for  Data  Set  One. 
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Figure  5.3.  Continued. 


the  factors  from  the  EBP  since  the  negative-going  features  have  been  almost  completely 
eliminated  from  the  spectral  and  time-mode  factors.  The  time  mode  of  Factor  Four  (water 
Raman  scatter)  exhibits  a  slight  bimodal  character. 

Three-mode  NNALS  (3M-NNALS)  was  also  executed  starting  with  EBP  factors  and 
with  the  same  convergence  criterion.  The  series  converged  after  93  iterations,  a  minor 
reduction  compared  with  3M-ALS.  The  results  of  this  calculation  are  in  Figure  5.4.  Not 
surprisingly,  the  3M-NNALS  factors  are  practically  indistinguishable  from  the  factors 
presented  in  Figure  5.3. 

Fluorescence  lifetimes  were  computed  with  the  phase  plane  method  of  Demas^^  from 
the  time-mode  profiles  shown  in  Figure  5.4.  The  Raman  scatter  from  WTM  1  was  used  as 
an  estimate  of  instrument  response  function,  £(r) .  The  results  were  7  ns  for  fluorene,  40 

ns  for  naphthalene,  and  140  ns  for  pyrene,  which  compared  favorably  with  independent 
determinations  of  7  ns,  36  ns,  and  127  ns,  respectively,  for  these  species  in  air-saturated 
water.  The  fluorescence  spectra  also  correlated  well  with  the  pure  component  spectra  at 
this  resolution. 

The  rank  5  decomposition  provided  an  intriguing  aspect  to  the  analysis.  A  contaminant 
from  an  unknown  source  appeared  as  a  fifth  factor.  The  results  of  the  EBP  and  the 
3M-NNALS  for  the  additional  factor  are  displayed  side-by-side  in  Figure  5.5,  It  is  not  clear 
at  this  time  what  the  source  of  the  fifth  factor  is  since  it  only  appears  in  one  WTM,  at 
approximately  440  nm  in  WTM  5  in  Figure  5.1.  The  lifetime  of  this  contaminant  is  very 
short,  approximately  2  ns. 
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Figure  5.4.  Continued. 
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Figure  5.5.  Rank  5  decomposition  of  Data  Set  One  -  Fifth  factor. 


The  3M-ALS  sequence  required  854  iterations  to  converge  for  the  rank  5  resolution 
when  starting  with  the  factors  from  the  EBP,  whereas  the  3M-NNALS  required  only  1 16 
iterations  to  converge,  starting  from  the  same  factors. 

The  other  resolutions,  ranks  three  and  six,  all  produced  factors  recognizable  as  the  three 
analytes  after  3M-ALS  or  3M-NNALS.  The  rank  three  resolution  coalesced  the  fluorene 
and  Raman  factors,  which  are  the  most  similar  of  the  four  chemical  components  in 
wavelength  and  in  fluorescence  decay  time.  In  the  TLD  for  assumed  rank  six,  all  of  the 
factors  from  the  rank  five  TLD  were  obtained.  In  addition,  a  small  magnitude, 
structureless,  short  lifetime  noise  factor  arose,  which  is  the  expected  result  when  a 
decomposition  of  rank  higher  than  the  actual  array  rank  is  sought. 

5.1.1.1.3  Comments  on  TLD  of  Data  Set  One 

When  confronted  with  relatively  simple  data,  such  as  that  contained  in  Data  Set  One, 
the  EBP  does  a  very  good  job  of  extracting  the  profiles  of  the  major  components  in  all  three 
modes.  3M-ALS  is  still  needed  to  make  further  refinements  in  the  profiles,  evidenced  by 
the  improvement  in  the  fluorene  factor  following  3M-ALS.  The  EBP’s  ability  to  extract 
weaker  components  is  not  as  good  as  3M-ALS.  One  may  even  have  overlooked  the 
mysterious  fifth  component  if  3M-ALS  had  not  been  performed. 

3M-NNALS,  on  the  other  hand,  did  not  represent  a  significant  improvement  over  the 
3M-ALS  result  with  respect  to  the  quality  of  the  factors.  Nonetheless,  3M-NNALS  is 
beneficial  because  it  provides  a  faster  route,  sometimes  substantially  faster,  to  convergence. 
It  is  not  clear  why  the  rank  five  3M-ALS  solution  took  so  long  to  converge,  but  the 
3M-NNALS  solution  apparently  avoided  that  problem. 


64 


3M-NNALS  is  undoubtedly  more  efficient  than  3M-ALS  at  avoiding  physically 
meaningless  regions  of  hyperspace.  3M-ALS  allows  the  solution  to  wander  through 
hyperspace  in  search  of  a  minimum.  We  have  observed  that  quite  often  two  corresponding 
profiles  for  a  factor  (time  and  wavelength,  say)  will  be  inverted.  It  is  possible  that  the 
phenomenon  of  "swamps"  described  by  Mitchell  and  Burdick^^  can  be  avoided  by 
3M-NNALS;  surely,  it  precludes  two  factor  degeneracies. 

Rank  estimation  is  another  area  requiring  comment.  Considering  the  Raman  scatter  and 
the  contaminants  as  independent  factors,  then  Data  Set  One  constitutes  a  rank  five  3-array. 
The  best  estimates  given  for  the  rank  seem  to  be  from  the  factor  indicator  function  and 
Malinowski’s  F-test  performed  on  the  summed  WTMs.  Other  than  the  estimate  for 
Solution  1,  there  were  no  rank  five  estimates  for  the  WTMs.  Even  Solution  5,  where  the 
contaminant  was  most  evident,  the  rank  estimate  was  “low.”  Summation  provides  some 
smoothing  or  signal  averaging  of  the  data,  which  may  be  enough  to  allow  recognition  of 
weaker  factors,  particularly  those  associated  with  processes  such  as  solvent  Raman 
scattering,  which  is  represented  in  all  the  WTMs  of  the  data  array. 

5.1.1.2  Data  Set  Two;  ATREEM 

Data  Set  Two  was  measured  on  a  single  aqueous  solution  that  contained  fluorene, 
phenanthrene,  naphthalene,  and  carbazole.  The  solution  was  prepared  by  adding  small 
amounts  of  stock  solutions  of  the  individual  components  in  methanol  to  3  mL  of  water  in  a 
quartz  cuvette  until  a  suitable  signal  level  was  reached.  A  TREEM  consisting  of  four 
WTMs  at  excitation  wavelengths  of  287  nm,  290  nm,  295  nm,  and  300  run  was  measured. 
The  data  interval  and  range  were  2  ns  over  250  ns  along  the  time  mode  and  2.5  nm  from 
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305  nm  to  405  nm  along  the  wavelength  mode.  Other  experimental  conditions  were 
identical  to  those  used  for  the  acquisition  of  Data  Set  One.  The  WTMs  that  comprise  the 
TREEM  for  Data  Set  Two  are  presented  in  Figure  5.6. 

Figure  5.6  reveals  that  there  is  much  less  variation  from  one  WTM  to  another  in 
comparison  to  Data  Set  One.  Note  that  the  variation  in  these  WTMs  arises  solely  from 
differences  in  the  spectral  response  due  to  changes  in  the  excitation  wavelength.  The 
concentrations  of  the  fluorescing  species  do  not  change  from  one  WTM  to  the  next.  One 
might  expect  the  estimation  of  the  rank  of  this  3-array  to  be  more  difficult  than  that  of  Data 
Set  One,  given  the  similarity  of  the  WTMs. 

5.1.1.2.1  Rank  Estimation  of  Data  Set  Two 

Rank  estimation  for  Data  Set  Two  was  performed  in  an  identical  maimer  as  Data  Set 
One.  The  results  of  these  analyses  are  presented  in  Table  5.3. 

Based  on  a  priori  knowledge  of  the  solution,  the  best  estimates  of  rank  for  this  data  are 
given  by  the  wavelength-mode  autocorrelation  and  Malinowski’s  F-test.  The  factor 
indicator  function  was  not  a  reliable  estimate  of  rank  estimation  for  any  of  the  WTMs  in 
this  data  set. 


Table  5.3.  Rank  estimates  for  Data  Set  Two. 
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^  Rank  is  equal  to  the  number  of  emitting  species. 

^  Malinowski's  factor  indicator. 

^  Malinowski's  F-test  for  significance  above  the  5%  level. 
^  Autocorrelation  coefficient,  values  above  +0.5. 
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Raman  scatter  was  not  considered  as  a  contributor  to  these  data.  Because  the  excitation 
wavelength  was  varied,  the  water  Raman  scattering  would  appear  at  a  different  position  in 
each  WTM.  If  the  Raman  scatter  were  included,  the  actual  rank  would  be  higher  by  a  value 
of  four.  However,  since  the  concentrations  of  the  fluorophores  was  elevated  to  optimize 
the  signal,  the  Raman  scatter  was  diminished  substantially  by  competition  for  incoming 
photons.  The  intensity  of  the  Raman  scatter  decreases  as  the  concentrations  of  the 
luminescent  species  are  increased.  This  behavior  can  be  observed  in  Figure  5.4,  Factor 
Four,  Solution-mode. 

5.1.1.2.2  TLD  of  Data  Set  Two 

EBP  decomposition  was  performed  on  Data  Set  Two  in  the  same  manner  as  for  Data 
Set  One.  Decompositions  were  performed  for  assumed  ranks  two  through  six.  Complex 
(imaginary)  factors  resulted  in  one  set  of  the  rank  three  factors  and  in  both  sets  of  the  rank 
four,  five,  and  six  factors  (Note:  this  EBP  produces  two  sets  of  factors,  one  corresponding 
to  eigenanalysis  on  the  G  matrices,  the  other  corresponding  to  eigenanalysis  on  the 
transpose  of  the  G  matrices).  Complex  factors  can  arise  during  eigenanalysis  on 
nonsymmetric  square  matrices  (see  Appendix  C),  which  is  generally  the  case  in  the 
approaches  required  here  (see  equations  C.12  and  C.13).  Li  et  al?^  devised  an  algorithm  to 
eliminate  complex  eigenvalues  during  EBPs  based  on  similarity  transforms. 

The  results  of  the  rank  four  EBP  utilizing  the  similarity  transform  method  of  Li  et 
are  displayed  in  Figure  5.7.  None  of  these  factors  can  clearly  be  assigned  to  any  one 
species.  They  appear  to  be  linear  combinations  of  the  individual  component  spectral  and 
time  profiles,  that  is,  they  are  rotated  together.  Each  of  the  various  rank  resolutions 
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Figure  5.7.  EBP  determined  factors  for  Data  Set  Two. 
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Factor  Four 


resulted  in  factors  similar  to  the  rank  four  resolution.  3M-ALS  is  obviously  required  in  this 


case. 

The  rank  four  3M-ALS  decomposition  of  Data  Set  Two  is  exhibited  in  Figure  5.8.  The 
respective  factors  from  Figure  5.7  were  used  as  starting  vectors  for  3M-ALS.  Convergence 
occurred  after  3846  iterations.  The  improvement  in  the  quality  and  appearance  of  the 
factors  is  noteworthy. 

Based  on  the  emission  wavelength-mode  and  time-mode  factors.  Factors  One  through 
Four  in  Figure  5.8  can  assigned  as  fluorene,  carbazole,  naphthalene,  and  phenanthrene, 
respectively.  Fluorescence  lifetimes  are  6  ns,  1 1  ns,  30  ns,  and  32  ns  for  Factors  One 
through  Four,  respectively.  They  compare  favorably  with  7  ns  for  fluorene  in  water,  15  ns 
for  carbazole  in  ethanol, 36  ns  for  naphthalene  in  water and  33  ns  for  phenanthrene  in 
water.^°  Lifetimes  were  calculated  with  the  phase  plane  method  with  a  separate 
measurement  of  scattered  laser  excitation  to  provide  E{t) . 

The  emission  wavelength-modes  for  Factor  Three  and  Factor  Four  indicate  that  the  two 
are  rotated  together.  Note  that  the  short  wavelength  portion  of  Factor  Four  has  a  very 
similar  structure  to  that  in  Factor  Three.  The  similarity  of  the  lifetimes  and  the  effects  of 
experimental  noise  are  each  likely  contributors  to  the  vector  rotation.  Even  if  the  3M-ALS 
process  is  allowed  to  continue  until  stricter  termination  criteria  are  met,  the  result  is 
essentially  the  same.  The  need  for  constraints  during  the  3M-ALS  sequence  is  more 
evident  for  Data  Set  Two  than  for  Data  Set  One. 

3M-NNALS  was  performed  on  the  data  using  the  factors  in  Figure  5.7  as  starting 
points.  The  sequence  converged  quickly,  taking  only  337  iterations.  Figure  5.9  contains 
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Figure  5.8.  3M-ALS  determined  factors  for  Data  Set  Two. 
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Figure  5.9.  Continued. 


the  rather  impressive  results.  Our  assignments  for  Factors  One  through  Four  are 
naphthalene,  phenanthrene,  carbazole,  and  fluorene,  respectively;  the  computed  lifetimes 
are  34  ns,  29  ns,  10  ns,  and  5  ns,  respectively.  Note  that  the  factors  have  not  been 
renumbered  from  those  given  by  the  algorithms  to  better  match  them  up  with  real  factors. 

The  3M-NNALS  factors  are  more  realistic  than  the  EBP  factors  and  are  also  an 
improvement  over  the  3M-ALS  result  since  the  rotation  of  the  naphthalene  and 
phenanthrene  factors  is  eliminated.  The  emission  mode  3M-NNALS  factors  are  compared 
with  fluorescence  spectra  measured  on  a  spectrofluorimeter  (Spec  2T2  Fluorolog)  in  Figure 
5.10.  The  agreement  is  excellent,  particularly  since  the  data  are  not  corrected  for  detector 
response  in  either  case. 

The  excitation  wavelength  modes  have  received  little  attention  thus  far.  This  is  because 
there  were  so  few  measurements  in  this  mode;  also  the  WTMs  were  not  corrected  for  laser 
power  variation  as  a  function  of  excitation  wavelength.  This  mode  provided  variations  in 
fluorescence  intensity  for  the  luminescent  species,  although  we  acknowledge  the  potential 
of  excitation  wavelength  as  an  additional  identifying  agent. 

Resolutions  for  the  other  assumed  ranks  were  less  satisfactory  than  the  rank  four 
results.  The  rank  three  resolution  produced  factors  similar  to  Factors  One,  Three,  and  Four 
in  Figure  5.8.  The  residual  hypersurface  of  the  rank  three  3M-NNALS  decomposition 
apparently  contained  MLO.  Some  starting  vectors  led  to  two  reasonable  factors  while 
others  led  to  only  one.  This  is  not  too  surprising  since  the  rank  three  resolution  is  rank 
deficient. 
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Intensity  (no  units) 


Factor  One  -  Naphthalene 


Factor  Two  -  Phenanthrene 


Factor  Three  -  Carbazole  Factor  Four  -  Ruorene 


Emission  Wavelength  (nm)  Emission  Wavelength  (nm) 


Figure  5.10.  Comparison  of  resolved  factors  to  pure  species  spectra.  Pure  component 
spectra  are  the  crossed  solid  lines.  Emission  wavelength-mode  factors  are  solid  hnes. 
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The  3M-NNALS  rank  five  decomposition  produced  reasonable  factors  for  fluorene, 
phenanthrene,  and  carbazole.  The  naphthalene  factor  appeared  to  be  rotated  with  the 
phenanthrene  factor.  However,  the  decomposition  did  not  generate  the  expected  noise 
factor.  Instead,  a  fifth  factor  emerged  with  a  lifetime  and  spectrum  similar  to  carbazole. 

The  rank  six  decompositions  failed  to  converge  after  more  than  eight  hours  of 
computational  time  for  both  3M-ALS  and  3M-NNALS. 

5.1.1.2.3  Comments  on  Data  Set  Two 

Data  Set  Two  represents  a  real  challenge  for  three-mode  decomposition,  especially  for 
EBP  and  3M-ALS.  The  rank  four  decomposition  by  3M-NNALS  is  highly  encouraging. 
Meaningful  profiles  were  recovered  in  the  spectral  and  time  domains,  and  the  number  of 
iterations  required  to  reach  convergence  was  reduced  over  the  unconstrained  procedures. 

In  all  probability,  the  most  significant  impediment  to  a  good  3M-ALS  solution  is  the 
presence  of  real  factors  with  very  similar  fluorescence  lifetimes.  The  prevention  of  factor 
rotations,  such  as  those  which  occurred  between  the  naphthalene  and  phenanthrene  factors, 
underlies  the  potential  of  3M-NNALS  and  similar  analyses,  such  as  restricted  Tucker 
models.^^ 

The  complex  solutions  to  the  EBP  analyses  point  to  what  may  be  an  instrumental 
challenge  in  the  acquisition  of  this  type  of  time-domain  data.  Sanchez  and  Kowalski 
suggested  that  complex  solutions  in  the  EBP  results  are  a  consequence  of  deviations  from 
trilinearity.  We  have  observed  that  the  dye  laser  pulse  duration  can  vary  with  wavelength. 
If  it  does,  then  the  fluorescence  decay  profiles  of  components  whose  lifetimes  are  not  much 
longer  than  the  laser  pulse  duration  are  similarly  affected,  and  the  WTMs  are  distorted.  One 
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way  to  overcome  this  problem  in  the  future  would  be  to  deconvolve  the  data  before 
three-mode  decomposition,  for  example,  by  the  exponential  series  method  for  fitting 
multiexponential  decays.^"*  Alternatively,  one  could  attempt  modifications  to  the  dye  laser 
itself  to  make  the  pulse  durations  as  constant  as  possible. 

Predictions  of  the  ultimate  potential  for  rank  estimation  from  TREEMs  is  difficult  with 
the  limited  data  presented  here.  A  priori  information  and  the  quality  of  the  rank  four 
3M-NNALS  decomposition  suggest  that  the  data  are  truly  rank  four.  However,  the 
emergence  of  an  apparent  fifth  factor  in  the  rank  five  decomposition  and  the  abundance  of 
complex  EBP  solutions  hint  at  a  rank  ambiguity.  The  estimate  of  rank  five  from  the  F-test 
may  be  further  support  for  such  a  conclusion.  Even  so,  the  F-test  operating  on  the  summed 
WTMs  still  provides  a  reasonable  place  to  start. 

It  is  unclear  why  the  rank  six  decomposition  failed  to  converge.  This  behavior  is 
unusual,  especially  for  the  3M-NNALS,  and  we  cannot  offer  an  explanation  at  this  time. 
5.1.2  Global  Analysis  of  Three-mode  Data 

The  power  of  global  analysis  lies  in  its  ability  to  simultaneously  model  some  physical 
observable  across  a  number  of  independent  measured  data  sets.  Thus,  it  provides  excellent 
opportunities  for  signal  averaging  and  can  describe  a  large  set  of  data  with  relatively  few 
parameters. 

5.1.2.1  Global  Analysis  of  Data  Set  One 

Data  Set  One  was  introduced  in  Section  5. 1.1.1.  There  is  a  particular  aspect  of  this  data 
set  which  allows  it  to  be  easily  analyzed  using  global  analysis.  The  Raman  scatter  in  the 
first  WTM  provides  an  excellent  system  response  function,  E(t) ,  for  convolution.  Since 
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global  analysis  essentially  performs  curve  fitting  using  the  convolution  of  E{t)  and  the 
exponential  decay  of  the  fluorophore,  y(r) ,  collecting  a  high  quality  E{t)  is  crucial. 

Global  analysis  was  conducted  using  the  Levenberg-Marquardt  search  algorithm  which 
is  packaged  in  the  Optimization  Toolbox  for  Matlab*.  Convolutions  were  performed  by 
numerical  integration  of  equation  2.1 1  using  a  polynomial  integration  program  written  in 
this  laboratory.  After  each  set  of  lifetime  parameters  was  estimated,  a  least  squares  fit  of 
the  data  was  performed  to  extract  the  solution-mode  and  wavelength-mode  parameters. 

The  3-array  was  unfolded  as  an  iX(yjt)  matrix  (see  Appendix  B)  with  the  i-mode 

representing  time,  the  y-mode  representing  solution,  and  the  ^-mode  representing  emission 
wavelength.  Termination  was  based  on  reaching  a  minimum  in  the  norm  of  the  residual 
matrix. 

To  model  these  data  based  on  five  components,  global  analysis  was  conducted  using  a 
zero  lifetime  and  four  nonzero  lifetimes.  The  excitation  profile  was  used  without 
convolution  as  a  zero  lifetime  component  to  accommodate  scatter.  Unfortunately,  a  global 
analysis  solution  to  a  five-component  problem  could  not  be  reached.  Attempts  generally 
resulted  in  two  very  similar  lifetimes,  usually  at  the  extremes.  Consequently,  solution  of 
the  least  squares  problem  leading  to  the  solution-wavelength-mode  was  poorly  defined  (i.e., 

if  A  represents  the  time-mode  matrix,  (A'A)“*  was  nearly  singular),  and  it  would  have  the 

properties  similar  to  a  two-factor  degeneracy  in  3M-ALS.  A  four-factor  global  analysis 
solution  was  pursued,  in  light  of  the  failure  of  the  five-factor  model. 

The  results  of  a  four-factor  (a  zero  lifetime  and  three  nonzero  lifetimes)  global  analysis 
are  exhibited  in  Figures  5.1 1  through  5.14.  Components  represented  by  the  factors  are  in 
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Figure  5.1 1.  Factor  One  global  analysis  result  for  Data  Set  One.  The  upper  plot  is  the  least 
squares  estimate  of  the  solution-wavelength-mode  factors.  The  lower  plot  is  the  instrument 
response  function. 
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Figure  5. 12,  Factor  Two  global  analysis  result  for  Data  Set  One.  The  upper  plot  is  the  least 
squares  estimate  of  the  solution-wavelength-mode  factors.  The  lower  plot  is  the 
convolution  of  the  instrument  response  with  the  decay  for  a  6.2  ns  lifetime. 
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Figure  5.13.  Factor  Three  global  analysis  result  for  Data  Set  One.  The  upper  plot  is  the 
least  squares  estimate  of  the  solution-wavelength-mode  factors.  The  lower  plot  is  the 
convolution  of  the  instrument  response  with  the  decay  for  a  36.4  ns  lifetime. 
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Figure  5. 14.  Factor  Four  global  analysis  result  for  Data  Set  One.  The  upper  plot  is  the 
least  squares  estimate  of  the  solution-wavelength-mode  factors.  The  lower  plot  is  the 
convolution  of  the  instrument  response  with  the  decay  for  a  133  ns  lifetime 


order:  Raman  scatter,  fluorene,  naphthalene,  and  pyrene.  The  lower  portion  of  the  figures 
contains  the  time-mode  part,  and  the  upper  portion  contains  the  least  squares  estimate  of 
the  solution-wavelength-mode.  In  Figure  5.1 1,  in  the  region  of  Solution  Five,  one  can  see 
the  contribution  of  the  contaminant.  It  can  be  seen  to  a  lesser  degree  in  the  same  region  in 
Figure  5.12.  Taken  together,  the  two  signals  can  be  fit  to  a  lifetime  of  2.8  ns,  almost 
identical  to  the  3M-ALS  result.  Why  this  factor  could  not  be  extracted  is  not  clear. 
Perhaps  3M-ALS  and  3M-NNALS  are  more  robust  techniques  for  extracting  factors  that 
make  relatively  small  contributions  to  the  data. 

The  results  for  the  remaining  factors  are  quite  good.  The  lifetimes  for  the  fluorene, 
naphthalene,  and  pyrene,  7  ns,  36  ns,  and  133  ns,  respectively,  are  in  good  agreement  with 
other  measurements  performed  in  this  laboratory  (see  Section  5. 1.1. 1.2).  They  differ  from 
those  obtained  from  3M-ALS  because  global  analysis  enforces  the  model  in  equation  2.1 1 
during  the  fitting  process,  while  3M-ALS  fits  the  data  to  the  trilinear  model  without  regard 
to  lifetime. 

Another  noteworthy  point  about  these  fits  is  the  behavior  of  the  scatter  or  zero-lifetime 
factor.  As  one  moves  from  Solution  One  to  Solution  Seven  in  the  top  portion  of  Figure 
5.1 1,  there  is  an  increasing  negative  contribution  of  the  zero-lifetime  factor.  This  effect 
becomes  markedly  more  prominent  when  a  discrete  convolution  is  performed  in  place  of  a 
numerical  integration.  This  effect  has  been  observed  in  this  laboratory  on  several 
occasions. 

Commonly,  the  problem  is  treated  as  a  time  shift  phenomenon  and  is  easily 
accommodated  in  this  way.  However,  it  has  also  been  observed  in  this  laboratory  that  a 
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small  time  shift  has  the  same  effect  as  a  short  lifetime  contribution.  One  possible  cause  of 
this  behavior  is  the  numerical  treatment  of  the  data;  another  is  that  the  response  of  PMT 
may  be  nonlinear  on  the  rising  edge. 

Carefully  comparing  the  time-mode  profiles  in  Figure  5.4  with  those  in  Figures  5.1 1, 
5.13,  and  5.14,  a  small  difference  in  the  rising  edges  can  be  detected.  The  global  analysis 
profiles  consistently  rise  faster  than  the  3M-NNALS  profiles.  This  may  be  a  contributing 
factor  to  the  problem  of  the  negative  signal  contribution;  and  it  may  have  the  same  cause. 
All  of  these  matters  are  under  investigation. 

5.1.2.2  Global  Analysis  of  Data  Set  Two 

Efforts  to  perform  global  analysis  on  Data  Set  Two  using  the  techniques  described  in 
Section  5. 1 .2. 1  were  unsuccessful.  Three  was  the  maximum  number  of  lifetimes  which 
could  be  determined.  Attempts  to  obtain  four  lifetimes  consistently  failed  due  to 
difficulties  identical  to  those  encountered  in  the  five-lifetime  analysis  of  Data  Set  One. 

The  three  lifetime  global  analysis  did  not  produce  clearly  recognizable  factors.  It  also 
was  extremely  sensitive  to  starting  point.  This  result  was  not  too  surprising,  given  the 
results  of  the  rank  three  decomposition  from  3M-ALS  and  3M-NNALS. 

Considering  the  success  of  the  rank  four  3M-NNALS  solution  for  this  data,  we  wrote 
and  performed  a  three-mode  global  analysis  algorithm.  The  algorithm  searched  lifetime 
parameter  space  in  a  nonlinear  fashion,  however,  for  each  new  lifetime,  the 
wavelength-mode  and  solution-mode  were  solved  by  a  two-mode  NNALS  routine,  rather 
than  solving  them  together. 
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The  four-lifetime  three-mode  global  analysis  produced  excellent  results.  Both  spectral 
and  temporal  profiles  were  nearly  indistinguishable  from  the  results  in  Figure  5.9. 

Lifetimes  extracted  were  5  ns,  1 1  ns,  35  ns,  and  38  ns,  for  the  fluorene,  carbazole, 
phenanthrene,  and  naphthalene  factors,  respectively. 

This  new  method  of  global  analysis  has  produced  encouraging  results.  Its  major 
drawback  is  that  it  requires  extensive  iteration. 

5.2  ProHIe  Extraction  of  Fuel  Fluorescence 

An  evaluation  of  the  performance  of  a  prototype  of  the  ROST™  system  in  comparison 
to  a  nitrogen  laser-based  system  designed  for  environmental  analysis  was  performed  as  part 
of  the  DoD  Site  Characterization  and  Analysis  Penetrometer  System  (SCAPS)  program. 
The  testing  was  conducted  at  the  Naval  Command,  Control  Ocean  Surveillance  Center, 
Research,  Development,  Test  and  Evaluation  Division  (NRaD)  during  February,  1994. 
During  the  evaluation,  numerous  data  sets  were  acquired  to  assess  the  capability  of  the  two 
fluorescence  systems  to  detect  and  identify  fuels  on  sod  matrices  and  to  quantitate 
contamination  levels.  The  following  four  fuels  were  selected  for  the  test  on  the  basis  of 
their  widespread  distribution  at  military  sites  throughout  the  country:  diesel  fuel-marine 
(DFM),  summer  grade  diesel  fuel  (DF),  unleaded  gasoline  (UG),  and  JP-4  jet  aircraft  fuel 
(JP4). 

In  one  part  of  the  evaluation,  participants  were  challenged  to  test  the  ability  of  the 
systems  to  identify  unknown  fuels.  Hence,  four  fuels  were  placed  onto  three  different  soil 
matrices  in  different  concentrations.  Sod  matrices  also  were  chosen  for  their  diversity;  the 
first  was  sand,  nonabsorbent  with  a  light  background;  the  second  was  from  Columbus  AFB, 


87 


MS,  (CAFB)  and  was  highly  absorbent  and  nonreflective;  the  third  was  from  China  Lake 
NAS,  CA,  (CLNAS)  and  was  also  absorbent  and  nonreflective.  Individual  fuel-soil 
mixtures  were  loaded  into  sets  of  soil  sample  holders  fitted  with  a  sapphire  window.  The 
loaded  sample  holders  rotated  during  measurement  to  minimize  photolysis  while  spectra 
were  simultaneously  measured  with  the  ROST™  prototype  and  the  nitrogen  laser  systems. 
Spectral  responses  of  the  12  combinations  of  the  four  fuels  on  the  three  soils  at  a 
concentration  of  3000  ppm  served  as  a  calibration  set.  As  a  blind  test  of  the  capabilities  of 
the  two  systems,  24  unknowns  were  prepared  in  concentrations  of  1000  or  10,000  ppm  and 
measured. 

The  fluorescence  data  collected  by  the  ROST™  prototype  were  acquired  for  290  nm 
excitation.  Twenty-one  emission  wavelengths  evenly  spaced  between  300  nm  and  500  nm 
were  monitored.  Decays  were  measured  at  1  ns  intervals  over  a  120  ns  range.  Fiber-optic 
light  delivery  and  collection  occurred  over  a  50  m  long  silica-clad  silica  two-fiber  probe. 
The  core  diameter  was  365  |j,m.  A  representative  WTM  for  each  of  the  fuel  products  is 
displayed  in  Figure  5.15. 

The  data  were  corrected  in  the  time-domain  for  the  wavelength  dependence  of  the 
transit  time  of  light  in  silica  optical  fiber,  t/o,  using^^ 


Vo  ~  „  ’ 
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where  I  is  the  optical  fiber  length,  and  v  is  the  group  velocity  of  light  in  silica  given  by; 
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Figure  5.15.  WTMs  of  four  fuels  on  soil.  Contours  represent  levels  of  equal  fluorescence 
intensity.  The  weakest  intensity  is  represented  by  the  outermost  contour  in  each  WTM. 
Each  WTM  here  has  15  contours,  however,  not  all  have  the  same  intensity  vdues.  This 
effectively  scales  the  WTMs  for  display  purposes. 


where  c  is  the  speed  of  light  in  vacuo,  and  n(A)  is  the  index  of  refraction  for  silica.  The 


index  of  refraction  was  determined  with  the  dispersion  equation:^® 


.  c,A^  C3A"  C5A"  I2 

iH - ! - 1 - 2 - 1 - 5 _ 

where  A  is  wavelength,  and  the  constants,  c,-,  are 


n(A)  = 


(5.3) 


Cl  =  0.6961663, 

C2  =  0.0684043, 

C3  =0.4079426, 

C4  =0.1162414, 

C5  =  0.8974794, 
and  Cg  =  9.86161. 

Data  points  are  lost  at  the  lower  and  upper  time  limits  due  to  the  correction.  There  were  107 
time  increments  after  correction.  Figure  5.16  contains  the  same  data  as  in  Figure  5.15  after 
the  correction  for  group  velocity  time  shift.  Note  that  the  data  in  Figure  5.16  have  a  more 
“square”  appearance  than  Figure  5.15. 

After  the  evaluation,  it  was  noted  from  the  soil  background  measurements  that  several 
of  the  sapphire  windows  contained  fluorescent  impurities  (vide  infra),  a  conclusion  thaf 
was  verified  by  direct  measurements  on  the  windows  themselves.  A  WTM  contour  plot  for 
one  of  the  fluorescent  sapphire  windows  is  displayed  in  Figure  5.17.  Unfortunately,  no 
record  was  kept  of  which  of  the  18  different  sapphire  windows  was  used  for  any  particular 
measurements  made  during  the  week-long  exercise. 

The  data  were  prepared  for  multimode  analysis  by  arranging  the  12  WTMs  of  the 
known  fuel-soil  mixtures,  the  24  WTMs  of  the  unknown  fuel-soil  mixtures,  and  four 
background  WTMs  (three  with  and  one  without  window  fluorescence)  into  a  3-array  of 
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Figure  5.16.  Time  shift  corrected  WTMs  of  four  fuels  on  soil.  Contours  represent  levels  of 
equal  fluorescence  intensity.  The  weakest  intensity  is  represented  by  the  outermost  contour 
in  each  WTM.  Each  WTM  here  has  15  contours,  however,  not  all  have  the  same  intensity 
values.  This  effectively  scales  the  WTMs  for  display  purposes. 
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dimension  21  by  107  by  40.  An  alternative  analysis  scheme,  which  wUl  be  discussed  in 
Section  5.2.3,  reduced  the  data  order  to  a  2-array  by  integration  along  the  time-mode. 
Integration  was  accomplished  simply  by  summing  along  the  time-mode  of  the  WTMs. 
Since  the  time  intervals  are  evenly  spaced,  this  procedure  is  identical  to  a  trapezoidal  rule 
integration  and  results  in  a  21  by  40  matrix  of  time-integrated  spectra. 


Wavelength  (nm) 

Figure  5.17.  Background  fluorescence  including  contaminated  sapphire  window. 


5.2.1  Rank  Estimation  of  Fuel  Data 

Rank  estimation  was  carried  out  on  the  fuel  fluorescence  data  in  the  same  manner  as 
applied  to  the  two  previous  data  sets.  Rank  was  also  estimated  from  the  matrix  of 
time-integrated  spectra.  Analyses  were  performed  on  the  data  both  with  and  without  the 
four  background  WTMs  to  evaluate  their  effect.  They  are  presented  in  Table  5.4;  results  on 
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the  individual  WTMs  have  been  omitted  for  brevity,  and  the  time-mode  autocorrelation 
estimate  was  not  pursued. 


Table  5.4.  Rank  estimates  for  fuel 
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Indicates  whether  only  36  fuel  WTMs  or  4  background  WTMs  were  included  ' 
Malinowski’s  factor  indicator. 

^  Malinowski's  F-test  for  significance  above  the  5%  level. 

Autocorrelation  coefficient,  values  above  +0.5. 

Based  on  the  results  for  the  two  previous  data  sets,  it  is  tempting  to  concentrate  on  the 
rank  estimates  given  by  the  F-test  and  the  wavelength-mode  autocoirelation.  These 
indicate  a  rank  of  around  four  for  this  data. 

The  factor  indicator  function  appears  to  be  unstable  with  respect  to  the  form  of  the  data. 
This  IS  probably  a  result  of  the  change  in  the  dimension  of  the  data  matrix,  since  the 
dimensions  are  used  in  the  computation  of  the  indicator  values.  The  summed  WTMs 
matrix  is  21  by  107,  while  the  time-integrated  matrix  is  21  by  40  or  21  by  36.  The  trend  for 
the  factor  indicator  function  is  to  increase  as  the  larger  dimension  is  decreased.  This  may 

beaweakness  of  the  factorindicator  function.  As  Mahnowski^  pointed  out,  the  factor 

indicator  function  is  an  empirical  function  and  should  be  used  with  caution. 


In  light  of  the  fact  that  the  fuels  are  known  to  be  complex  mixtures  with  many  different 
PAHs,  the  results  presented  in  Table  5.4  are  surprising.  It  is  possible  that  there  is  only  a 
small  number  of  emitter  groups  and  the  individual  compounds  within  the  groups  have  very 
similar  emission  spectra  and  lifetimes.  Other  considerations  include  the  level  of 
experimental  noise,  species  interactions,  matrix  effects,  and  group  velocity  correction  errors 
for  the  fiber  optics.  All  of  these  factors  may  contribute  to  the  rank  estimate  being  far  lower 
than  expected  for  these  data. 

5.2.2  TLD  of  Fuel  Fluorescence 

Mitchell  and  Burdick^^  have  proposed  classifying  mixtures  via  trilinear  decomposition 
of  three-mode  data.  For  example,  consider  a  3-array  consisting  of  WTMs  of  many  different 
fuel  samples.  Two  of  the  modes  from  the  trilinear  decomposition  would  represent 
wavelength  and  fluorescence  decay  time,  and  the  third  mode  would  then  represent  a 
mixture  or  fuel-type  mode.  The  relative  proportions  of  the  factors  in  the  class-mode  could 
be  used  as  descriptors  in  a  classification  scheme,  thereby  exploiting  the  second-order 
advantage.  One  could  not  only  classify  the  fuel,  but  also  identify  its  major  components. 

To  obtain  a  unique  decomposition  for  3-mode  data,  the  matrices  which  form  the 
underlying  triple  product  of  the  data  must  meet  the  condition  of  Equation  B.40,  namely  that 
the  sum  of  the  k-ranks  of  the  matrices  must  be  greater  than  or  equal  to  two  times  the  rank  of 
the  3-array  plus  two.  Thus,  if  even  one  of  those  matrices  has  two  identical  columns,  a 
unique  decomposition  is  not  possible  with  3M-ALS. 

Is  it  likely  that  the  fuel  data  collected  have  two  identical  or  proportional  columns  in  the 
factor  matrices  that  make  up  the  3-arTay?  Since  the  fuels  are  complex  mixtures  containing 
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tens  or  even  hundreds  of  fluorescent  species  and  since  there  are  only  four  fuels,  it  is  likely 
that  the  mixture-mode  factor  matrix  has  k-rank  of  one.  Therefore,  we  expect  that  the 
3-array  cannot  be  uniquely  decomposed  with  3M-ALS.  Even  considering  the  matrix  effects 
of  the  different  soils,  which  might  contribute  to  variability  among  the  fluorescence  from 
different  species,  the  number  of  fluorophores  is  likely  too  large  to  obtain  a  unique 
decomposition. 

Another  consideration  in  analyzing  fuel  emission  is  the  possibility  of  nonlinear 
spectroscopic  effects  arising  from  interactions  between  species.  Some  of  the  potential 
sources  of  deviations  from  trilinearity  for  individual  fluorophores  are  quenching,  energy 
transfer,  and  reabsorption  of  fluorescence. 

The  more  factors  that  are  sought,  the  longer  is  the  required  computational  time. 

Trilinear  decomposition  of  this  data  set  by  3M-ALS  or  3M-NNALS  took  as  long  as  a  day 
(using  a  120  MHz  Pentium®  computer  and  employing  strict  convergence  criteria)  for 
assumed  rank  of  ten  or  greater,  which  is  impractically  long  for  many  applications.  Given 
these  complications,  what  can  be  expected  or  gained  from  a  trilinear  decomposition  here? 

Eastwood  has  suggested  that  the  sources  of  fluorescence  in  complex  mixtures,  such  as 
fuels,  could  be  grouped  into  classes  of  fluorophores.  The  fluorescence  spectral  and  lifetime 
data  offer  some  encouragement.  For  example,  single  ring  compounds  typically  emit  in  the 
short  wavelength  range,  say  280  nm  to  300  nm,  and  their  lifetimes  are  also  usually  short  in 
air-saturated  aqueous  solution.  Naphthalene  and  its  derivatives  have  intermediate 
wavelength  emission  (300  nm  -  360  nm)  and  lifetimes  (ca.  25-35  ns)  emission.  And  the 
higher  molecular  weight  PAHs  tend  to  emit  at  even  longer  wavelengths;  however,  the 
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trends  in  lifetime  are  less  obvious  for  the  PAHs.  Perhaps  this  is  the  best  that  can  be 
expected  from  fluorescence  data  of  fuel  products  in  view  of  the  results  from  the  simple 
solutions.  Recall  that  in  the  four-component  case  of  Data  Set  Two,  we  had  to  employ 
3M-NNALS  to  obtain  a  satisfactory  result. 

If  the  goal  of  the  decomposition  is  extraction  of  the  fluorescence  factors  of  all  species 
contributing  to  the  signal,  trilinear  decomposition  is  probably  not  a  realistic  approach;  and 
more  traditional  analytical  methods  such  as  GC-MS  or  HPLC  should  be  considered. 
However,  if  the  goal  is  to  find  factors  which  adequately  describe  the  fuels  while  taking 
advantage  of  the  third-order  nature  of  the  data  and  generate  good  descriptors  in  the  process, 
then  3M-NNALS  may  prove  useful. 

Given  the  previous  rationalizations,  we  attempted  to  keep  the  rank  of  the  decomposition 
as  low  as  possible,  but  large  enough  to  model  the  data  and  the  sapphire  window 
interference.  Resolutions  of  ranks  three  through  eight  were  obtained  using  EBP 
decomposition.  All  resulted  in  complex  (imaginary)  factors.  The  similarity  transform 
method  was  employed  to  eliminate  the  complex  factors.  The  resulting  factors  had  generally 
large  negative  amplitude  contributions  with  little,  if  any,  physical  significance.  A  factor 
clearly  associated  with  the  window  fluorescence  did  not  emerge  from  the  analysis. 

The  EBP  computations  were  followed  by  3M-ALS.  The  results  from  the  unrestricted 
3M-ALS  were  also  not  very  encouraging.  The  calculations  tended  to  produce  two  or  three 
very  similar  factors.  Although  these  did  not  always  meet  the  definition  of  a  two-factor 
degeneracy,  they  were  usually  highly  correlated  in  all  three  modes  and  were  physically 
meaningless.  A  factor  indicative  of  the  sapphire  window  fluorescence  could  be  identified 
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in  some  resolutions,  but  not  in  others.  Computation  times  were  very  long,  as  greater  than 
50,000  iterations  were  required  to  reach  convergence  in  several  cases. 

3M-NNALS  was  performed  for  ranks  three  through  seven  using  the  results  of  the 
similarity  transform  method  decomposition  as  starting  vectors.  The  results  from  these 
calculations  appeared  more  like  real  factors.  However,  rather  than  having  smooth 
asymptotic  transitions  to  the  baseline,  the  mode  features  often  tended  to  get  zeroed  out  or 
“cut  off’  as  they  reached  zero  intensity  or  concentration  for  many  factors.  This  is  an 
unfortunate  artifact  of  the  non-negativity  algorithm.  On  the  other  hand,  the  fact  that  none 
of  the  decompositions  performed  with  this  method  required  over  300  iterations  for 
convergence  is  encouraging. 

The  rank  five  resolution  by  3M-NNALS  that  produced  the  smallest  residual  is 
presented  in  Figure  5.18.  Unfortunately,  3M-NNALS  (at  least  in  its  current 
implementation)  is  very  susceptible  to  trapping  in  local  minimum  solutions.  At  least  two 
different  resolutions  were  determined  for  the  rank  five  case.  Both  resolutions  produced 
reasonable  representations  of  the  window  contribution,  and  the  other  factors  were  similar. 

Three  background  samples  which  obviously  contained  window  fluorescence  were 
analyzed  separately,  also.  The  background  WTMs  with  window  fluorescence  were  trivial 
to  analyze  and  proved  to  be  rank  one,  based  on  the  appearance  of  rank  one  through  three 
3M-NNALS  resolutions.  The  factor  for  the  interference  is  displayed  for  comparison  in  the 
last  frame  of  Figure  5.18. 

The  rank  five  resolution  was  chosen  from  among  all  of  the  decompositions  because  it 
provided  the  best  match  to  the  window  factors  obtained  from  the  3M-NNALS  analysis  of 
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Figure  5.18.  Rank  five  3M-NNALS  decomposition  factors  for  fuel  fluorescence 
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Figure  5.18.  Continued. 
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Figure  5.18.  Continued. 


just  the  background  WTMs  themselves.  Clearly,  the  fifth  factor  of  the  rank  five  resolution 
represents  the  interference.  It  highly  correlates  with  the  background  resolution  in  both  the 
wavelength  and  time  modes;  the  angle  cosines  of  the  overlap  are  0.9973  and  0.9984  in  the 
respective  modes.  The  results  for  the  higher  rank  resolutions  are  similar  to  those  of  rank 
five.  The  rank  four  resolution  did  not  provide  a  clear  identification  of  the  window  factor. 

It  is  unfortunate  that  the  features  appear  to  get  “cut  off’  for  a  number  of  the  factors. 
Nevertheless,  the  window  fluorescence  can  be  identified  and  removed  from  each  sample 
WTM  by  subtracting  the  triple  product  of  the  vectors  for  the  three  modes  of  Factor  Five  in 
Figure  5. 18  from  the  3-array. 

Sample-mode  Factors  One  through  Four  were  used  as  descriptors  in  a  classification 
scheme  (see  Section  5.2.3). 

5.2.3  Analysis  of  Fuel  Fluorescence  With  Self-modeling  Curve  Resolution 

A  fundamental  question  about  the  multimode  data  generated  by  a  ROST™  or  similar 
system  is  whether  all  this  information  is  really  needed.  Certainly  first-order  data,  in  the 
form  of  spectra,  can  be  collected  very  efficiently  using  a  charge  coupled  device  (CCD)  or 
photo  diode  array  (PDA)  detector.  The  time-integrated  data  created  during  rank  estimation 
represents  such  first-order  data.  Such  data  can  be  utilized  to  compare  the  merits  of  higher 
order  data  with  lower  order  data. 

The  time-integrated  data  were  subjected  to  a  self-modeling  curve  resolution  (SMCR) 
scheme  using  a  two-mode  NNALS  algorithm.  Decompositions  were  obtained  for  ranks 
three  through  seven.  The  wavelength-mode  results,  in  many  cases,  were  similar  to  those  of 
the  wavelength-mode  from  3M-NNALS.  However,  we  found  that  the  extracted  factors 
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were  very  sensitive  to  starting  points.  An  important  aspect  of  the  analysis  was  extracting 
the  sapphire  window  component.  In  many  cases  of  profile  extraction  here,  a  single  window 
factor  did  not  emerge  from  the  calculation.  This  was  the  case  for  our  rank  five  minimum 
residual  result. 

5.2.4  Comments  on  TLD  of  Fuel  WTMs 

3M-NNALS  is  the  only  feasible  way  to  analyze  WTMs  from  complex  mixtures  such  as 
fuels.  EBP  and  3M-ALS  were  not  able  to  obtain  physically  meaningful  factors  from  these 
data,  and  the  latter  was  far  too  time  intensive. 

Assignments  of  the  factors  to  groups  of  compounds,  such  as  naphthalenes  and 
substituted  naphthalenes,  etc.,  is  premature;  and  further  analysis  of  fuel  fluorescence  in 
needed.  Performing  3M-NNALS  on  simulated  fuel  mixtures  which  are  chemically  and 
spectroscopically  similar  to  fuels,  but  of  known  composition,  may  aid  in  this  process.  Rank 
estimation  performance  should  also  benefit  from  this  type  of  research. 

The  true  signifrcance  of  the  TLD  technique  is  its  ability  to  identify  and  remove 
interferences  from  the  data  and  generate  descriptors  of  the  fuels.  3M-NNALS  was  able  to 
isolate  a  factor  closely  resembling  the  sapphire  window  interference  in  the  rank  five 
resolution,  an  accomplishment  that  two-mode  analysis  could  not  match.  How  successfully 
and  completely  the  technique  removed  the  window  factor  cannot  be  characterized  for  these 
data  because  the  true  amount  of  window  interference  is  not  known.  But  in  the  analysis  that 
follows,  we  can  speculate  on  its  performance  based  on  classification  results  with  and 
without  the  window  component.  More  work  needs  to  be  performed  in  this  area  also. 
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5.3  Classification  of  Fuels  Using  Fluorescence  Data 

One  element  of  the  field  laser  performance  evaluation  was  to  determine  the  ability  of 
the  instruments  to  classify  the  24  unknown  fuel-soil  mixtures  using  the  12  known  fuel-soil 
mixtures  as  a  training  set. 

To  separate  the  fuel  groups,  we  used  Fisher’s  discriminant  function,  described  in 
Section  4.1.2,  to  generate  linear  discriminants  with  which  to  classify  the  unknowns. 
Multivariate  normal  distributions  and  equal  covariance  matrices  for  the  four  fuels  were 
assumed.  Prior  probabilities  were  equal  under  the  conditions  of  the  evaluation.  Although 
the  classification  algorithm  presented  in  Section  4.1.2  is  reasonably  straightforward,  the 
challenge  is  to  find  the  most  accurate  descriptors  of  the  fuels  within  the  data. 

In  addition  to  using  the  mixture-mode  factors,  descriptors  for  the  classification 
algorithm  were  obtained  in  several  other  ways.  One  method  used  the  mixture-mode 
elements  from  three-mode.  Another  method  used  the  mixture-mode  PFA  factors 
(computed  using  S  VD)  of  the  time-integrated  spectra  as  descriptors.  Each  of  the  methods 
used  was  evaluated  on  its  ability  to  predict  fuel  type  of  the  24  unknowns.  Normalization  of 
the  descriptor  vectors  was  performed  in  all  cases  to  minimize  the  effects  of  concentration 
and  the  soil  matrix.  The  normalization  consisted  of  scaling  the  vector  comprising  the 
descriptors  of  each  fuel  to  unit  length.  However,  the  stage  of  the  classification  process  at 
which  normalization  was  performed  depended  on  how  the  descriptors  were  obtained. 

5.3.1  3M-NNALS  Mixture-mode  Factors  as  Descriptors 

The  elements  of  a  3M-NNALS  mixture-mode  factor  represent  the  contribution  of  the 
corresponding  wavelength-time  dyad  to  the  WTMs.  Thus,  the  third  element  of 


103 


mixture-mode  Factor  One  in  Figure  5.18,  for  example,  is  the  contribution  of  the 
wavelength-time  dyad  of  Factor  One  needed  to  model  the  WTM  of  the  third  mixture 
(which  happens  to  be  a  DFM  on  sand  combination).  Each  element  of  each  mixture-mode 
factor  is  a  descriptor  of  one  of  the  fuels  (or  window).  A  rank  five  decomposition  yields 
five  descriptors  for  each  fuel;  a  rank  six  decomposition,  six;  and  so  forth. 

Table  5.5  contains  the  classification  results  utilizing  all  of  the  factors,  including  the 
window  factors,  from  3M-NNALS  decomposition  for  ranks  four  through  seven. 

To  eliminate  the  window  fluorescence  from  the  data,  the  obvious  window  factor  was 
identified  for  each  rank  decomposition,  and  its  mixture-mode  factor  was  omitted  as  a 
descriptor.  The  descriptor  vectors  were  normalized  and  used  in  classification.  The 
classification  results  with  the  window  contributions  omitted  are  in  Table  5.6. 

The  classification  rates  were  good  for  all  rank  decompositions,  with  and  without  the 
window  factors.  In  all  cases  of  misclassification,  the  mixtures  corresponded  to  the  lowest 
fuel  concentration  (1000  ppm)  in  the  evaluation  samples.  In  fact,  the  most  consistently 
misclassified  mixture.  Unknown  M,  was  a  low  concentration  of  JP4,  the  weakest  emitter  of 
the  fuels,  and  was  mixed  on  the  CAFB  soil,  which  was  the  poorest  soil  matrix  in  terms  of 
fluorescence  response.  Support  for  this  last  statement  comes  from  the  fact  that  10  of  the  12 
misclassifications  in  Table  5.5  occurred  for  CAFB  soil  samples;  in  Table  5.6,  eight  of  the 
10  misclassifications  occurred  on  the  CAFB  soil.  Unknown  M  had  a  very  weak  emission 
and  the  poorest  S/N  of  any  sample  in  the  study. 

We  suspect  that  the  second  most  misclassified  mixture.  Unknown  T,  has  a  very  large 
sapphire  window  interference  in  proportion  to  its  total  fuel  emission.  This  is  also  the  case 
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Table  5.5.  Classification  of  unknown  fuels  using  3M-NNALS  factors 
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Table  5.6.  Classification  of  unfaiown  fuels  using  3M-NNALS  factors  omitting  suspected 
window  factor  _ 
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for  Unknown  E,  although  to  a  lesser  degree. 

The  integrated  fluorescence  intensities  versus  unknown  mixture  are  displayed  in  Figure 
5.19.  The  estimated  window  fluorescence  has  been  removed  from  the  fluorescence  before 
integration  using  the  rank  five  3M-NNALS  window  factor.  Included  in  Figure  5.19  are  the 
integrated  window  fluorescence  estimates  for  each  mixture.  Note  the  intensities  and 
window  components  of  the  misclassified  mixtures.  It  is  obvious  here  which  unknowns 
suffer  from  the  greatest  interference. 

In  general,  fuel  fluorescence  intensity  follows  the  trend:  DFM  is  slightly  stronger  than 
DF,  which  is  stronger  than  UG,  which  is  much  stronger  than  JP4.  The  trend  for  soil  effects 
based  on  fuel  fluorescence  intensity  on  the  respective  soil  is  fuel-sand  intensity  is  much 
greater  than  fuel-CLNAS  soil  mixture  intensity,  which  is  greater  than  fuel-CAFB  soil 
intensity. 

It  is  surprising  and  disappointing  that  the  number  of  misclassifications  does  not 
improve  more  when  the  suspected  window  component  is  omitted  from  the  set  of 
descriptors.  A  disturbing  aspect  of  this  analysis  is  that  three  classifications.  Unknowns  K 
and  N  in  the  Four  Factor  colunrn  and  Unknown  T  in  the  Five  Factors  column,  went  from 
correct  with  the  window  descriptor  included  to  incorrect  when  it  was  omitted.  It  is  also 
difficult  to  establish  a  pattern  for  misclassification  in  the  presence  of  window  fluorescence. 
Window  fluorescence  may  make  UG  appear  as  DFM,  as  in  the  case  of  Unknown  R,  or  the 
reverse,  as  in  the  case  of  Unknown  E. 

Apparently,  the  two  conditions,  weak  emission  and  large  proportional  window 
fluorescence  contribution  to  total  signal,  are  the  most  significant  causes  of 
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misclassification.  However,  it  is  not  clear  yet  that  the  results  of  Table  5.6  are  the  best 
which  can  be  achieved. 

Conceivably,  the  trilinear  model  employed  here  has  failed  to  model  subtleties  in  the 
data  which  are  important  in  classification.  This  may  be  a  fault  with  using  the  3M-NNALS 
factors  as  descriptors,  since  they  are  not  the  true  least  squares  estimates  of  the  data.  It  is 
possible  that  PFA,  with  its  orthonormal  basis  set  spanning  the  true  factor  space,  may  be 
able  to  generate  better  descriptors  for  the  weaker  intensity  fuel  mixtures. 


w 
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Figure  5.19.  Unknown  fuel  fluorescence  and  sapphire  window  fluorescence.  Integrated 
fluorescence  intensities  for  unknown  mixtures  are  represented  by  o’s.  Integrated  window 
fluorescence  intensities  are  represented  by  +  signs. 
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5.3.2  Three-mode  PFA  Mixture-mode  Factors  as  Descriptors 

The  rank  five  sapphire  window  factor,  Factor  Five  of  Figure  5.18,  was  subtracted  from 
the  WTMs  of  the  12  training  mixtures  and  24  unknowns.  The  resulting  3-array  was 
analyzed  using  the  TUCKERl  (eigenanalysis  on  the  unfolded  3-array  without  performing 
ALS),  three-mode  PFA. 

The  mixture-mode  eigenvectors  were  multiplied  by  their  corresponding  eigenvalues  and 
used  as  descriptors.  Before  discrimination  and  classification,  the  descriptors  were 

normalized  as  in  the  Section  5.3.1.  The  results  for  various  numbers  of  PFs  are  presented  in 
Table  5.7. 

The  results  of  the  classification  here  are  excellent.  One  interpretation  that  can  be  maHp. 
is  that  the  factors  from  the  3M-NNALS  do  not  adequately  represent  the  data.  This  is 
probably  a  result  of  substantial  deviations  from  trilinearity  of  the  fuel  WTMs.  However, 
the  trilinear  model  seems  to  be  adequate  to  fit  and  extract  the  sapphire  window 
fluorescence. 

The  TUCKALS3  algorithm  was  also  employed  to  generate  descriptors  with  identical 
results. 

5.3.3  Principal  Factors  of  Time-mtegrated  Spectra  as  Descriptors 

The  WTMs  of  the  12  training  mixtures  and  the  24  unknowns  were  time  integrated  both 
before  and  after  correction  for  group  velocity  shift  of  the  fiber  optics.  The  data  were 
arranged  into  a  36  by  21  matrix  and  subjected  to  PFA  using  SVD  such  that  the  left  singular 
vectors  spanned  the  row  space  (i.e.,  the  mixture-mode  space)  of  the  matrix  and  the  right 
singular  vectors  spanned  the  column  space  (i.e.,  the  wavelength-mode  space)  of  the  matrix. 
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Table  5.7.  Classification  of  unknown  fuels  using  three-mode  PFA  mixture-mode  PFs 


MMIliwIii:! 
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DF 

DF 

DF 

DF 
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JP4 
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Total 

24 

24 

24' 

23 

Correct 

100% 

100% 

100% 

96% 

- - - 1 - ; - 1 - 1 - 1 - 

Used  rank  five  3M-NNALS  decomposition  to  remove  sapphire  window  (interference) 


factor. 
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Descriptors  were  generated  by  taking  the  product  of  the  left  singular  vector  matrix  with 
the  singular  value  (diagonal)  matrix.  Then,  they  were  entered  into  the  linear  discrimination 
and  classification  algorithm.  The  capsule  results  of  the  analyses  are  contained  in  Table  5.8. 


Table  5.8.  Number  of  correct  fuel  classifications  with  PFA  of  time-inteprated  spectra 


Source^:;!! 

PFs 

Raw 

WTMs 

20 

21 

22 

21 

Corrected 

WTMs 

20 

21 

23 

23 

‘  Spectra  generated  from  summing  (time-integrating)  WTMs  before  correcting  for  group 
velocity  shift  of  fiber  optics  (Raw  WTMs)  and  after  correction  (Corrected  WTMs). 


The  misclassified  fuels  from  these  analyses  were  all  low  concentration.  The  same  fuels 
were  misclassified  in  the  three  PF  and  four  PF  cases.  This  was  not  the  case  with  the  five 
PF  and  six  PF  cases.  Table  5.8  is  not  intended  to  illustrate  a  problem  with  time  integration 
of  spectra  or  the  use  of  PFA  on  spectra  or  the  group  velocity  shift  Correction  of  the  data. 
Indeed,  the  spectra  in  the  two  cases  differ  only  minutely  due  to  losses  of  information  on  the 
edges  of  the  time-mode.  The  PFs  also  have  only  minor  differences  even  up  to  rank  six. 

Table  5.8  points  out  a  possible  problem  with  the  classification  algorithm.  When  the 
sample  size  of  the  standards  (training  set)  is  small  compared  to  the  number  of  descriptors, 
the  covariance  matrix  estimates,  S,  of  equation  4.8,  can  become  highly  variable  (and 
unstable  upon  inversion)  and  result  in  poor  classification.  The  number  of  standards  per 
class  here  is  three.  When  the  number  of  descriptors  exceeds  three,  the  classification  may  be 
suffering.  The  results  in  the  PF  5  and  PF  6  columns  of  Table  5,8  may  be  a  manifestation  of 
this  phenomenon. 
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The  discrepancy  illustrated  in  Table  5.8  might  have  another  cause.  Recall  that  the  rank 
estimation  for  the  matrix  of  integrated  spectra  was  four  or  five.  Eigenvectors  greater  than 
four  or  five  may  be  noise  eigenvectors  and  disrupt  the  classification.  Although  the  variance 
accommodated  by  the  higher  eigenvalues  is  very  small,  their  effect  may  be  more  prominent 
on  the  weaker  fluorescing  mixtures. 

After  removing  the  window  factors  as  in  Section  5.3.2,  the  WTMs  were  time  integrated 
and  factor  analyzed  with  SVD.  The  significant  factors  were  then  used  in  the  discrimination 
and  classification  algorithms  to  classify  the  unknowns.  Table  5.9  contains  the 
classifications  based  on  the  3M-NNALS,  rank  five  window  removal. 

The  results  in  Table  5.9  are  superb.  Apparently  choosing  to  integrate  and  using 
orthogonal  factors  in  the  classification  scheme  was  an  excellent  choice.  Since  these  results 
are  perfect,  the  inference  can  be  made  that  classification  of  these  fuels  does  not  require 
three-mode  data.  However,  the  analyses  required  to  get  to  this  point  speak  for  themselves 
concerning  the  three-mode  nature  of  the  data. 

Subtraction  of  3M-NNALS  window  factors  other  than  the  rank  five  decomposition 
were  also  performed.  The  classification  results  using  PFs  from  the  time-integrated  spectra 
in  these  cases  were  not  as  good  as  the  rank  five  window  removal  case.  However,  they  were 
as  good  or  better  than  the  3M-NNALS  mixture-mode  results  presented  in  Table  5.6.  Again, 
deviations  from  trilinearity  in  the  data  give  advantages  to  the  orthogonal  factors  to 
represent  the  data.  Clearly  though,  the  choice  of  window-mode  factor  is  very  important; 
thus,  so  is  the  choice  of  number  of  factors  in  the  3M-NNALS  model. 


112 


Table  5.9.  Classification  of  unknown  fuels  using  PFA  of  time-integrated  spectra 
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Used  rank  five  3M-NNALS  decomposition  to  obtain  sapphire  window  (interference) 
factor. 


Analogous  classification  was  performed  on  the  two-mode  data  after  subtraction  of  the 
window  factors  generated  with  SMCR.  The  results  were  an  improvement  over  those  in 
Table  5.8,  but,  in  general,  not  as  strong  as  in  Table  5.9.  Subtraction  of  the  rank  four  SMCR 
window  factor  resulted  in  one  perfect  classification  using  four  PFs. 


113 


5.3.4  Comments  on  Fuel  Classification 


This  data  set  was  particularly  difficult  to  deal  with  because  of  the  interference. 
Certainly,  greater  care  could  have  been  exercised  in  the  selection  of  the  sapphire  windows 
used  for  the  measurement.  On  the  other  hand,  this  type  of  problem  is  probably  highly 
representative  of  those  that  will  be  encountered  during  in  situ  analysis,  whether  at  a 
hazardous  waste  site  or  in  human  tissue.  The  effort  expended  here  should  facilitate 
analyses  on  actual  unknowns. 

Correct  classification  of  the  fuels  was  possible  after  removing  the  interference  and 
modeling  the  data  with  orthogonal  factors  to  generate  descriptors.  Major  considerations  of 
classification  were  the  method  of  profile  extraction,  selection  of  window  factor,  choice  of 
descriptors,  and  the  limited  size  of  the  training  set. 

Based  on  these  results  and  the  results  presented  in  Section  5.2,  we  believe  that  the  best 
approach  to  profile  extraction  of  fuel  fluorescence  is  3M-NNALS.  While  3M-NNALS  is 
susceptible  to  MLO,  it  is  less  so  than  SMCR. 

Future  studies  should  be  conducted  to  include  larger  training  sets  with  more  fuel  types. 
They  should  include  various  levels  of  oxygen  to  examine  the  effects  of  fluorescence 
quenching  on  profile  extraction  and  classification. 

5.4  Cluster  Analysis 

Cluster  analysis  could  be  extremely  useful  in  distinguishing  different  types  of 
luminescent  materials,  such  as  fuels  or  tissue  types,  without  prior  knowledge  of  types  of 
materials  present.  Such  a  scheme  could  be  applied  to  a  contaminated  site  containing 
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multiple  fuels  or  to  a  collection  of  tissue  specimens  containing  healthy  and  diseased  tissues 
or  tissue  from  different  species. 

The  goal  of  this  section  is  to  evaluate  the  viability  of  HCA  to  perform  clustering  of 
fuels.  In  addition,  we  would  like  to  find  the  best  measure  of  similarity  to  be  employed  in 
the  HCA  algorithm. 

As  in  most  of  the  classification  schemes  used  in  this  work,  HCA  was  performed  on  data 
corrected  for  the  sapphire  window  contribution  by  subtracting  the  window  factor  from  the 
3M-NNALS  rank  five  decomposition.  The  HCA  algorithm  described  in  Section  4.2  was 
performed  with  Matlab*  using  a  program  written  in  this  laboratory.  Numerous  approaches 

were  taken  for  both  generating  of  the  similarity  (dissimilarity)  matrices  and  forming 
linkages  during  clustering. 

Clustering  was  performed  using  all  of  mixture-mode  PFs  employed  in  the 
classifications  of  Section  5.3  as  descriptors  of  the  fuel  data.  Both  “city  block”  and 
Euclidean  metrics  were  tested  as  dissimilarity  measures.  Single,  complete,  and  average 
linkage  methods  were  investigated  for  forming  linkages.  All  of  these  descriptors,  the 
3M-NNALS  factors,  PFs  from  time-integrated  spectra,  and  PFs  from  three-mode  PFA, 
performed  clustering  very  well. 

The  average  linkage  dendrogram  of  the  descriptor-vector-normalized  3M-NNALS 
factors,  excluding  the  sapphire  window  factor,  is  displayed  in  Figure  5.20.  The  Euclidean 
metric  was  used  to  generate  the  distance  matrix.  Using  a  distance  of  around  0.4  as  a  group 
separation  cutoff,  only  Unknowns  T  and  M  are  incorrectly  grouped.  In  fact,  they  do  not 
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Figure  5.20.  Average  linkage  dendrogram  of  fuels  using  3M-NNALS  factors.  Normalized 
mixture-mode  factors  from  rank  five  3M-NNALS  decomposition  of  fuels  with  sapphire 
window  (interference)  factor  omitted.  Similarity  measure  is  Euclidean  distance  between 
descriptor  vectors.  The  known  spectra  are  indicated  by  their  fuel  abbreviations  and  the 
unknowns  by  their  letter  designators. 
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group  with  any  of  the  major  fuel  groups.  This  is  not  too  siuprising  given  that  these  same 
two  unknowns  were  misclassified  using  these  descriptors  (see  Table  5.6). 

Figure  5.21  portrays  the  average  linkage  dendrogram  of  the  fuels  using  the  first  four  left 
singular  vectors  of  the  36  by  21  matrix  of  time-integrated  spectra  multiplied  by  their 
respective  eigenvectors.  These  are  the  same  data  that  were  used  in  the  4  PF  column  of 
Table  5.9.  With  a  distance  cutoff  of  0.2,  only  Unknown  M  is  not  grouped  with  any  major 
group.  Otherwise,  the  groupings  are  excellent. 

PFs  from  three-mode  PFA,  used  in  4  PF  column  of  Table  5.7,  were  used  to  produce  the 
average  linkage  dendrogram  in  Figure  5.22.  Here  the  groupings  are  more  ambiguous, 
although  still  very  good.  If  one  chooses  a  distance  cutoff  of  0.5,  one  of  the  JP4  standards 
and  Unknown  A  are  excluded  from  the  JP4  grouping.  Choosing  a  cutoff  of  1.0  lumps  the 
UG  and  DF  fuels  into  a  single  group.  Unknown  M  is  incorrectly  grouped  into  the  UG 
cluster  in  any  event. 

When  dealing  with  spectral  vectors,  the  angle  cosine  is  one  of  the  more  intuitive 
measures  of  similarity  because  it  indicates  the  proximity  of  two  vectors  in  space  and  is 
therefore  a  measure  of  correlation.  The  angle  cosine  calculation  has  normalization  built 
into  it;  only  the  direction  of  the  vector  is  considered,  and  the  magnitude  is  unimportant. 

For  these  reasons,  the  remainder  of  this  section  will  deal  with  clustering  in  which  the  angle 
cosine  (or  its  value  subtracted  from  unity)  is  taken  as  the  measure  of  similarity. 

Figure  5.23  contains  the  single  linkage  dendrogram  using  one  minus  the  angle  cosine  of 
the  spectra  as  the  measure  of  similarity.  There  is  excellent  separation  of  the  groups  based 
upon  these  criteria.  Only  one  mixture.  Unknown  M,  did  not  group  as  expected,  which  is 
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Figure  5.21.  Average  linkage  dendrogram  of  fuels  using  PFA  of  time-integrated  spectra. 
First  four  mixture-mode  PFs  of  normalized  spectra.  Sapphire  window  (interference)  factor 
removed  using  rank  five  3M-NNALS.  Similarity  measure  is  Euclidean  distance  between 
descriptor  vectors.  The  known  spectra  are  indicated  by  their  fuel  abbreviations  and  the 
unknowns  by  their  letter  designators. 
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Figure  5.22.  Average  linkage  dendrogram  of  fuels  using  three-mode  PFA  factors. 
Normalized  mixture-mode  factors.  Sapphire  window  (interference)  factor  removed  using 
rank  five  3M-NNALS.  Similarity  measure  is  Euclidean  distance  between  descriptor 
vectors.  The  known  spectra  are  indicated  by  their  fuel  abbreviations  and  the  unknowns  by 
their  letter  designators. 
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Figure  5.23.  Single  linkage  dendrogram  of  time-integrated  fuel  spectra.  The  angle  cosine 
of  the  spectral  vectors  subtracted  from  unity  was  used  as  the  measure  of  similarity. 
Sapphire  window  (interference)  factor,  obtained  from  rank  five  3M-NNALS 
decomposition,  was  removed  before  integration.  The  known  spectra  are  indicated  by  their 
fuel  abbreviations  and  the  unknowns  by  their  letter  designators. 
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not  surprising  based  on  the  results  in  the  preceding  section.  One  can  easily  see  four  groups 
using  a  similarity  index  of  approximately  0.015  as  a  cutoff. 

Based  on  the  similarity  index  values  in  Figure  5.23,  one  might  get  the  mistaken 
impression  that  the  spectra  are  almost  identical  except  for  very  subtle  variations.  This 
arises  because  the  single  linkage  method  finds  the  minimum  distance  or  similarity  between 
groups  as  they  are  established.  In  reality,  the  angle  cosines  (correlations)  between  different 
fuel  spectra  can  be  quite  small,  as  low  as  0.45  for  DFM  and  JP4.  Figure  5.23  does  correctly 
indicate  the  general  relative  similarity  between  the  fuel  spectra,  i.e.,  UG  and  DF  have  the 
most  similar  spectra,  and  DFM  is  more  similar  to  UG  and  DF  than  JP4. 

Figure  5.24  contains  the  complete  linkage  dendrogram  for  the  same  data  depicted  in 
Figure  5.23.  Note  that  the  groups  have  greater  separation  than  in  the  single  linkage  model. 
Unknown  M  also  does  not  appear  to  be  as  well  separated  from  the  other  groups.  These 
effects  result  from  the  choice  of  maximum  distance  or  similarity  between  groups  as 
opposed  to  the  minimum  in  the  single  linkage  method.  The  movement  of  Unknown  M  in 
relation  to  the  other  groups  as  the  method  is  changed  indicates  of  the  spread  of  the  angle 
cosines  within  the  groups. 

While  the  groups  seem  to  have  a  large  variance,  given  the  large  changes  in  the 
similarity  indices  from  single  linkage  to  complete  linkage  and  the  movement  of  Unknown 
M,  the  groups  are  well-defined  and  well-separated.  This  interpretation  arises  from  the 
formation  of  the  same  groups  regardless  of  method  chosen.  If  different  groups  were  formed 
as  the  method  of  HCA  changed,  which  is  common,  then  a  very  judicious  choice  of  groups 
would  be  required. 
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Figure  5.24.  Complete  linkage  dendrogram  of  time-integrated  fuel  spectra.  The  angle 
cosine  of  the  spectral  vectors  subtracted  from  unity  was  used  as  the  measure  of  similarity. 
Sapphire  window  (interference)  factor,  obtained  from  rank  five  3M-NNALS 
decomposition,  was  removed  before  integration.  The  known  spectra  are  indicated  by  their 
fuel  abbreviations  and  the  unknowns  by  their  letter  designators. 
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Figure  5.25  contains  the  average  linkage  dendrogram  for  the  data  depicted  in  Figures 
5.23  and  5.24.  Once  again,  the  groups  are  nicely  separated  and  distinct.  A  similarity  index 
cutoff  of  approximately  0.4  separates  the  four  fuels,  leaving  Unknown  M  isolated. 

Based  on  the  results  presented,  HCA  is  an  excellent  method  of  grouping  fuel  mixtures. 
The  average  linkage  method  of  clustering  is  probably  as  good  a  choice  as  any,  although  the 
single  and  complete  linkage  methods  should  be  examined  in  each  case.  The  angle  cosine  of 
the  time-integrated  fluorescence  spectra  is  an  excellent  measure  of  similarity,  and  it  is  easy 
to  apply  to  fluorescence  emission  data.  This  should  be  a  valuable  tool  in  site 
characterization. 

Of  course,  as  in  classification,  it  is  important  to  numerically  remove  any  interferences 
before  HCA.  If  the  sapphire  window  component  is  not  removed,  many  of  the  low 
concentration  mixtures  would  not  cluster  correctly. 
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Figure  5.25.  Average  linkage  dendrogram  of  time-integrated  fuel  spectra.  The  angle 
cosine  of  the  spectral  vectors  subtracted  from  unity  was  used  as  the  measure  of  similarity 
Sapphire  window  (interference)  factor,  obtained  from  rank  five  3M-NNALS 
decomposition,  was  removed  before  integration.  The  known  spectra  are  indicated  by  their 
fuel  abbreviations  and  the  unknowns  by  their  letter  designators. 
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6.  CONCLUSIONS  AND  RECOMMENDATIONS 


6.1  Summary  of  the  Present  Work 

This  research  represents  the  first  application  of  TLD  and  global  analysis  methods  to  a 
WTM-concentration  3-array  and  a  TREEM.  This  is  also  the  first  use  of  TLD  methods  on 
highly  complicated  data  in  the  form  of  fluorescence  of  fuels.  Previous  works  have  applied 
TLD  methods  to  EEM-concentration/quencher  level  arrays^*^’^^’^^  and  EEFAs^^"^^  of 
solutions  with  only  a  few  components.  TLD  of  WTM  based  data  gives  greater  flexibility 
than  EEM  based  data  since  it  can  be  fit  to  a  well-established  model,  exponential  decay, 
which  can  be  utilized  in  a  three-mode  global  analysis  scheme.  WTM  data  typically  have 
better  S/N  than  EEFA  data  because  optical  switches  used  in  phase-modulated  instruments 
do  not  provide  perfect  modulation. 

Global  analysis  of  fluorescence  decay  time-mode  data  has  been  used  to  decompose 
WTMs  of  two  component  mixtures;  however,  global  analysis  of  three-mode  data  has  not 
been  performed. 

The  efficacy  of  classification  and  clustering  algorithms  applied  to  reduced  forms  of 
three-mode  fuel  fluorescence  data  was  also  demonstrated. 

Analyses  of  Data  Set  One  illustrate  the  capabilities  of  various  rank  estimation  and 
profile  extraction  techniques  with  low  rank  data.  Data  Set  One  consists  of  seven  WTMs. 
The  WTMs  were  the  fluorescence  response  of  six  solutions  containing  varying  amounts  of 
fluorene,  pyrene,  and  naphthalene  in  water  and  a  water  blank.  Including  the  water  Raman 
scatter  as  a  component,  the  3-array  has  an  a  priori  rank  of  the  3-array  of  four. 
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Excellent  rank  estimates  of  the  matrix  summed  WTMs  of  Data  Set  One  were  obtained 
using  Malinowski’s  factor  indicator  function,  F-test,  and  first-lag  autocorrelation 
coefficient  of  the  wavelength-mode  eigenvectors.  The  first-lag  autocorrelation  coefficient 
of  the  time-mode  eigenvectors  did  not  give  good  rank  estimates.  However,  time-mode 
based  estimates  could  be  improved  by  using  a  cutoff  value  greater  than  +0.5,  which  was 
used  for  the  wavelength  mode.  This  is  the  first  use  of  matrix  summation  of  the  WTMs 
before  rank  estimation  reported.  It  represents  a  significant  improvement  of  the  capabilities 
of  rank  estimation  of  3-arrays. 

The  four  component  WTM-concentration  3-array  in  Data  Set  One  was  easily 
decomposed  with  an  EBP,  3M-ALS,  and  3M-NNALS.  3M-ALS  and  3M-NNALS  offered  a 
slight  improvement  to  the  EBP  result.  3M-NNALS  offered  a  minor  computation  time 
improvement  over  3M-ALS. 

Profile  extraction  using  global  analysis  on  Data  Set  One  provided  good  results  when 
seeking  three  nonzero  lifetimes  and  including  a  zero  lifetime  scatter  component.  However, 
global  analysis  did  not  perform  well  when  attempting  to  also  extract  the  profile  of  a  very 
weak  interference,  which  was  discovered  in  the  three-mode  analyses. 

Data  Set  Two  is  a  four-component  TREEM  of  a  single  solution  of  fluorene, 
naphthalene,  carbazole,  and  phenanthrene.  The  spectra  of  carbazole  and  phenanthrene  are 
heavily  overlapped,  and  the  lifetimes  of  naphthalene  and  phenanthrene  are  very  similar. 
These  factors  made  for  a  difficult  decomposition  problem. 

Performance  of  rank  estimation  methods  were  essentially  the  same  as  for  Data  Set  One. 
A  matrix  representing  the  compressed  3-array  was  produced  by  matrix  summation  of  the 
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four  WTMs  obtained  at  the  excitation  wavelengths  utilized  to  generate  the  TREEM. 
Excellent  rank  estimates  were  obtained  using  Malinowski’s  factor  indicator  function,  F- 
test,  and  first-lag  autocorrelation  coefficient  of  the  wavelength-mode  eigenvectors.  The 
first-lag  autocorrelation  coefficient  of  the  time-mode  eigenvectors  did  not  give  good  rank 
estimates  using  +0.5  as  a  cutoff  value. 

TLD  was  performed  on  this  data  using  an  EBP,  3M-ALS,  and  3M-NNALS.  The  EBP 
and  3M-ALS  were  unable  to  provide  realistic  factor  profiles  for  this  data.  This  was 
probably  the  result  of  the  similarity  in  lifetimes  of  naphthalene  and  phenanthrene. 
3M-NNALS  produced  very  good  results,  generating  recognizable  factors  in  one-tenth  of  the 
time  required  by  3M-ALS.  3M-NNALS  provides  substantially  improved  results  over 
3M-ALS,  a  notable  advancement  of  three-mode  decomposition. 

Global  analysis,  performed  using  the  methodology  employed  with  Data  Set  One,  failed 
to  produce  a  satisfactory  result  for  any  number  of  factors.  The  algorithm  either  failed  to 
obtain  recognizable  factors  or  terminated  with  two  nearly  identical  lifetimes.  The  latter 
condition  resulted  in  an  ill-conditioned  solution  of  the  excitation-emission  mode. 

A  hybrid  form  of  analysis,  which  combines  global  analysis  with  NNALS,  provided  very 
promising  results  that  were  nearly  identical  to  3M-NNALS.  However,  the  algorithm  was 
very  time  intensive.  This  utilization  of  the  time-mode  information  in  conjunction  with 
three-mode  analysis  is  unique.  We  believe  this  methodology  has  great  potential  as  an 
analysis  technique  because  there  is  a  model  in  the  time  mode  which  does  not  exist  in  the 
excitation  or  emission  modes. 
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The  final  set  of  data  examined  in  this  research  was  composed  of  WTMs  of  fuels  on  soil 
matrices.  This  data  contained  an  interference  in  the  form  of  sapphire  window  fluorescence. 
The  goal  of  analysis  was  to  use  TLD  to  remove  the  interference,  and,  in  the  process, 
generate  descriptors  of  the  data  that  could  be  used  in  a  classification  scheme. 

Rank  estimation  was  performed  for  this  data  in  the  same  manner  as  in  the  two  previous 
data  sets.  The  rank  estimates  were  smaller  than  expected,  considering  the  number  of 
fluorescing  species  in  fuels. 

TLD  using  an  EBP  was  not  successful,  yielding  complex  eigenvectors  or  chemically 
meaningless  factors.  3M-ALS  performed  better  than  the  EBP,  but  it  produced  multiple 
factor  degeneracies  and  took  a  great  deal  of  computation  time.  3M-NNALS  performed  best 
of  the  three  TLD  procedures.  It  generated  factors  which  may  be  interpreted  as  real  factors, 
although  no  assignment  was  attempted.  Computation  time  was  at  least  two  orders  of 
magnitude  faster  than  3M-ALS.  The  rank  five  decomposition  provided  the  best  match  for 
the  sapphire  window  factor. 

Classifications  of  the  unknown  fuels  using  the  mixture-mode  factors  from  3M-NNALS 
were  very  good.  Using  all  of  the  rank  five  factors  as  descriptors,  22  out  of  24  unknown 
fuels  were  correctly  classified.  Omitting  the  suspected  window  factor  and  using  the 
remaining  four  factors  as  descriptors,  22  out  of  24  unknowns  were  classified  correctly.  One 
of  the  unknowns  was  misclassified  in  both  cases.  In  all  cases  of  misclassification,  the  fuels 
were  low  concentration  and  were  usually  on  the  soil  matrix  with  the  poorest  fluorescence 
response. 
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Classification  rates  were  improved,  often  to  100%,  when  the  window  factor  was 
subtracted  fi'om  the  3-array.  After  performing  eigenanalysis  on  the  3-array  or  the 
time-integrated  spectra  of  the  fuels,  the  mixture-mode  eigenvectors  were  used  as 
descriptors  of  the  fuels. 

HCA  was  performed  utilizing  metrics  computed  firom  the  various  descriptors  described 
in  the  preceding  paragraph  and  one  minus  the  angle  cosine  as  measures  of  similarity.  HCA 
was  a  good  method  of  grouping  the  fuel  mixtures.  Using  the  angle  cosine  as  the  measure  of 
similarity  is  intuitive  and  gave  excellent  cluster  results.  Forming  linkages  with  the  all 
linkage  methods  worked  well,  and  each  should  be  utilized  for  each  application  of  HCA. 

6.2  Recommendations  for  Future  Work 

Future  efforts  in  this  area  should  be  directed  at  evaluating  more  fuels  and  fuel  types. 
These  analysis  techniques  may  be  able  to  detect  subtle  differences  in  fuels,  but  this  must  be 
explored. 

To  understand  the  complexities  of  fuel  component  interactions,  mixtures  of  simulated 
fuels  can  be  prepared  and  analyzed  using  these  fluorescence  techniques.  Simulated  fuels 
could  easily  be  prepared  in  a  progressive  fashion,  allowing  examination  of  the  interactions 
as  they  develop.  Using  this  methodology,  phenomena  such  as  the  low  rank  of  the  fuel 
WTMs  may  be  addressed. 

Oxygen  quenching  can  be  problematic  in  time-mode  fluorescence  spectroscopy.  In  the 
presence  of  oxygen,  the  intensity  of  the  fluorescence  is  diminished  and  the  lifetime 
shortened.  If  two  samples  of  the  same  fuel  undergo  different  levels  of  quenching,  their 
components  will  have  different  lifetimes,  but  identical  emission  modes.  This  will  make  for 
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a  difficult  analysis  problem  because  the  time-mode  will  have  a  higher  rank  than  the 
emission  mode.  Restricted  Tucker  models^’  are  designed  to  deal  with  problems  of  this 
type.  Current  programs  can  be  easily  modified  and  tested  to  perform  this  type  of  analysis. 

The  3M-NNALS  routine  works  by  setting  all  negative  elements  in  a  matrix  to  zero  after 
each  iteration.  This  rather  abrupt  method  may  lead  to  trapping  the  solution  in  a  local 
minimum.  A  more  amenable  procedure  would  involve  setting  the  negative  elements 
progressively  less  negative  each  iteration.  For  example,  reduce  the  negative  elements  by 
1%  on  the  first  iteration;  and  on  successive  iterations,  reduce  them  by  an  additional  1%.  In 
this  way,  the  unrealistic  solutions  are  guided  gently  back  into  the  region  of  real  solutions. 

Since  fluorescence  is  easily  capable  of  higher  orders  than  used  in  the  present  work, 
consideration  should  be  given  to  acquisition  and  analysis  of  fourth  order  data.  Excitation- 
emission-time-concentration  data  could  be  generated  at  the  present  time  with  no 
modifications  to  existing  equipment.  Advantages  of  this  type  of  data  over  third-order  data 
are  not  known.  However,  an  additional  mode  of  data  may  allow  degeneracies  in  one  mode 
to  be  lifted  in  another,  for  example,  quenching. 

Time-mode  data  is  unique  because  an  entire  vector  of  data  can  be  expressed  as  a  single 
variable,  lifetime,  for  a  single  component.  Greater  effort  should  be  directed  at  developing 
the  three-mode  global  analysis.  One  method  which  may  be  of  interest  would  force  the 
time-mode  to  fit  the  convolutional  model  of  equation  2.1 1  after  each  iteration  of  a  3M-ALS 
or  3M-NNALS  sequence.  This  may  generate  more  realistic  solutions  in  some  cases. 

These  physical  methods  and  fluorescence  analysis  techniques  should  be  applied  to  areas 
outside  of  environmental  analysis.  Medical  diagnostics  is  one  area  that  has  been  exploring 
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fluorescence.  Advances  in  environmental  analysis  may  be  readily  applied  to  diagnosis  of 
cancer  and  cardiovascular  disease  and  medical  monitoring. 

Classification  of  healthy  and  diseased  tissue  utilizing  third-  or  higher-order  data  with 
one  mode  of  data  as  a  healthy/diseased-mode  (analogous  to  the  fuel  type  mode)  is  easily 
conceivable.  This  type  of  analysis  would  facilitate  identification  of  the  biochemical 
indicators  of  disease. 
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APPENDIX  A.  DEFINITIONS  AND  NOTATION 


Throughout  this  thesis,  I  used  the  conventions  drawn  from  the  chemometric  literature 
for  displaying  vectors  and  matrices.^^’^*’^®'**  Scalar  quantities  are  represented  by  lower  case 
letters.  Scalars  which  are  elements  of  a  vector  or  matrix  are  labeled  by  subscript  indices 
indicating  row  and  column.  Column  vectors  are  indicated  by  bold  lower  case  letters,  e.g., 
a.  Row  vectors  are  denoted  by  a  prime,  e.g., a' .  Matrices  are  represented  by  bold  upper 
case  letters,  e.g.,  A.  Therefore, 

Uy  is  a  scalar  (the  subscript  denotes  the  {ij}th  element  of  the  matrix  A). 


is  the  /th  column  vector  (in  m-coordinate  space). 


a'  =  [fljj  0,2  •••  is  the  ith  row  vector  (in  n-coordinate  space). 


DeHnition  A.I.  (Inner  Product)  The  inner  product  (also  called  the  dot  product  or  scalar 
product)  of  two  vectors  is 


Definition  A.I.  Note  1.  The  inner  product  is  a  scalar,  and  it  requires  that  both  vectors 
have  the  same  number  of  elements. 
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DeHnition  A.l.  Note  2.  The  properties  of  the  inner  product  are  preserved  whether  the 
vectors  are  row  or  column  vectors,  i.e., 
ab  =  a'b  =ab' 

DeHnition  A.l.  Note  3.  (Vector  Norm)  The  length  (or  norm)  of  a  vector  is  given  by  the 
square  root  of  the  inner  product  of  the  vector  with  itself: 

IHI  =  =  ^  =  ^ 

Deflnition  A.2.  (Normal  Vector)  A  normal  vector  has  unit  length.  That  is,  for  the 
normal  vector  a, 

|a||  =  l- 

Deflnition  A.2.  Note  1.  A  vector  may  be  normalized  by  setting  its  length  equal  to  one. 
This  is  done  by  dividing  the  vector  by  its  length,  e.g.. 


a 


Deflnition  A.3.  (Orthogonal  Vectors)  Two  vectors  are  orthogonal  if  their  inner  product 
is  equal  to  zero  (0). 

Deflnition  A.4.  (Orthonormal  Vectors)  Two  vectors  are  orthonormal  if  their  inner 
product  is  equal  to  zero  and  their  lengths  are  equal  to  one  (1). 

Deflnition  A.5.  (Matrix  Dimension)  A  matrix  that  has  m  rows  with  n  elements  in  each 
row  is  said  to  have  matrix  dimensions  m  by  n.  A  convenient  shorthand  for  representing  the 
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m  by  n  matrix  A  is  ,  so 


A  —  ~ 

DeHnition  A.5.  Note  1.  In  many  texts,  this  property  of  a  matrix  is  called  the  order  of  the 
matrix.  However,  the  term  order  is  used  to  describe  a  variety  of  properties  of  matrices 
and  arrays  in  various  texts  and  publications  and  here. 

Definition  A.6.  (Matrix  Addition)  Matrix  addition  is  the  addition  of  the  corresponding 
elements  of  matrices  of  identical  dimensions.  For  example,  the  sum  of  matrices  A  3,^2 


Uji  a,2 

^21  ^22 

.^31  ^32. 


is  3  by  2  a  matrix. 


B 


3x2 


is 


i 

^12  “^^12 

^3x2"^®  3x2 

^21  "^^21 

^22  *^^22 

ajj  +i>3, 

^32  “^^32  . 

Definition  A.6.  Note  1.  Matrix  addition  of  matrices  of  different  dimensions  is 
undefined. 

Definition  A.6.  Note  2.  Vector  addition  is  completely  analogous  to  matrix  addition 
where  one  of  the  matrix  dimensions  is  one  (1). 

Definition  A.7.  (Scalar  Multiplication)  Multiplication  of  a  scalar  by  a  matrix  is  the 
product  of  each  element  of  the  matrix  by  the  scalar  multiplier.  The  product  of  a  and  B3,,2 

is 


‘aft,!  ab^2 

^®3x2 

^21  ^22 

z= 

ab^i  ab^ 

P^\  ^32 . 

_ab2i  abj2_ 
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Definition  A.7.  Note  1.  Scalar  multiplication  of  a  vector  is  completely  analogous  to 
scalar  multiplication  of  a  matrix. 

Definition  A.8.  (Matrix  Multiplication)  Multiplication  of  two  matrices,  denoted  AB,  is 
accomplished  by  computing  the  inner  product  of  the  ith  row  of  the  left  matrix  (A)with  the 
jth  column  of  the  right  matrix(B).  Matrix  multiplication  requires  that  the  two  matrices  be 
conformable,  that  is,  the  number  of  columns  of  the  left  matrix  must  equal  the  number  of 
rows  of  the  right  matrix.  Therefore,  if  A^^^  and  are  to  be  multiplied,  then  n  must 

equal  p.  If  n  ^  p ,  then  multiplication  is  not  defined.  For  example. 


AB  = 


‘2  4  3‘ 

2  0  f 

2  1  1 

zr 

1  4  3_ 

_2  3  1_ 

4+0+2  8+0+3  6+0+1 
2+8+6  4+4+9  3+4+3 


'6  11  7" 

16  17  10 

Definition  A.9.  (Matrix  Rank)  The  rank  of  a  matrix  is  the  maximum  number  of  linearly 
independent  column  (or  row)  vectors  in  the  matrix. 

DeHnition  A.9.  Note  1.  The  rank  of  a  matrix  is  the  true  underlying  order  of  that  matrix. 
Rank  is  also  defined  as  the  order  of  the  largest  matrix  with  a  nonzero  determinant 
formed  by  deleting  rows  and  columns  of  the  original  matrix.  The  maximum  rank  of  the 
matrix  A„,^  is  the  minimum  value  of  m  and  n.  The  rank  of  A„^  is  denoted  by 


(^mxn)‘ 

Definition  A.9.  Note  2.  A  matrix  whose  rank  is  less  than  minimum  of  its  dimensions  is 
said  to  be  rank  deficient. 

Definition  A.9.  Note  3.  A  zero  matrix  (a  matrix  all  of  whose  elements  are  zero)  has  rank 

0. 
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Definition  A.10.  (Square  Matrix)  A  square  matrix  has  the  same  number  of  rows  and 
columns.  Therefore,  the  matrix  A„,^„  is  square  if  m=n. 


Definition  A.ll.  (Matrix  Transpose)  The  transpose  of  a  matrix  is  obtained  by 
interchanging  the  rows  and  columns  of  the  matrix.  Transpose  is  indicated  by  a  prime,  so. 


If  A  = 


2  1  3‘ 
5  4  7_ 


then  A'  = 


'2  5* 
1  4 
_3  7_ 


Definition  A.12.  (Diagonal  Matrix)  A  diagonal  matrix  is  a  square  matrix  with  zero 
elements  except  along  the  principal  diagonal.  Any  matrix.  A,  is  diagonal  if  =0  for  all 


j ,  and  if  at  least  one  ^  0  when  i  =  j.  Therefore, 


2000 
0  300 
0000 
0006 


is  a  diagonal  matrix. 


Definition  A.13.  (Identity  Matrix)  An  identity  matrix,  I,  is  a  diagonal  matrix  with  ones 
along  the  principal  diagonal.  Since  the  identity  matrix  is  always  a  square  matrix,  its 
dimensions  are  indicated  by  a  single  subscript,  I„.  The  matrix  product  of  any  matrix  with 
the  identity  matrix  is  itself. 

'1  0---0' 

0  1  •••0 

1=  :  ;  .  : 

0  0-  1 
■^3x2l2~^3x2 
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Definition  A.14.  (Matrix  Inverse)  The  inverse  of  an  n  by  n  (i.e.,  a  square  matrix)  is 
another  matrix  such  that  the  product  of  the  two  matrices  is  an  identity  matrix.  The  inverse 
of  a  matrix  is  denoted  by  a  superscript  (-1).  An  invertible  matrix  is  said  to  be  nonsingular. 

A  A  —  T 

"3x3^  3x3  ~  ^3 

Definition  A.14.  Note  1.  A  rank  deficient  matrix  will  not  be  invertible  since  it  will  have 

a  determinant  equal  to  zero.  Rank  deficient  matrices  are  also  called  singular  matrices. 

Definition  A.15.  (Matrix  Pseudoinverse)  The  pseudoinverse  of  a  matrix  is  the  least 

squares  solution  to  a  set  of  linear  equations.  It  is  designated  by  a  superscripted  cross  (+). 

b  =  Ax 
A'b  =  A'Ax 

(A 'A)"’ AT)  =  (A'A)"‘ A'Ax 
(A 'A)"*  AT)  =  Ix  =  x 
A^b=x 

A"^  s(A'A)~'A's  pseudoinverse(A) 

Definition  A.15.  Note  1.  The  pseudoinverse  is  usually  used  when  A  is  a  rectangular 
matrix  (i.e.,  m  by  n),  since  A  has  no  true  inverse..  Also,  the  condition  m>n  must  be 
satisfied  or  the  matrix  (A 'A)  will  be  singular. 

Definition  A.16.  (Orthonormal  Matrix)  A  matrix  whose  column  vectors  are  mutually 
orthonormal  is  a  column- wise  orthonormal  matrix.  It  has  the  convenient  property: 

A'A=I. 

Likewise,  a  matrix  whose  row  vectors  are  mutually  orthonormal  is  a  row-wise  orthonormal 
matrix,  and 

AA'=I. 
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DeHnition  A.16.  Note  1.  (Orthogonal  Matrix)  A  matrix  whose  column  vectors  are 
mutually  orthogonal  is  a  colunm-wise  orthogonal  matrix.  It  has  the  property: 
A'A=D 


where  D  is  a  diagonal  matrix. 

Likewise,  a  matrix  whose  row  vectors  are  mutually  orthogonal  is  a  row-wise  orthogonal 
matrix,  and 
AA'=D. 

Definition  A.17.  (Matrix  Square  Root)  Let  Y  be  the  matrix  square  root  of  X,  then, 
YY=X. 

and, 

Y=XV^ 

Definition  A.18.  (Matrix  Norm)  The  norm  of  a  matrix  is  a  scalar  that  gives  a  measure 
of  the  magnitude  of  the  elements  of  the  matrix.  Here  the  matrix  norm  is  defined  as: 

W  =  (A.3) 


Definition  A.19.  (Angle  Cosine)  The  cosine  of  the  angle,  0 ,  between  two  vectors  a  and 


b  (also  called  the  angle  cosine)  is 

al) 


COS(t>  =  ■ 


IN 


(A.4) 


Definition  A.20.  (N-way  Array)  An  array  is  an  ordered  collection  of  scalars.  An  A^-way 
array,  or  N-array,  is  an  ordered  collection  of  scalars  which  have  indices  that  extend  info  N- 
dimensions.  A  matrix  is  a  2-array  and  is  a  rectangular  array  of  scalar  elements.  Likewise,  a 
3-array  is  a  rectangular  parallelepiped  array  (a  “box”  or  “cube”)  of  scalar  elements.  Higher 
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order  arrays  (“hyperboxes"  or  “hypercubes”)  extend  into  higher  dimensional  space  and  do 
not  have  a  geometrical  analog,  but  they  are  readily  conceivable  and  do  exist. 

Definition  A.20.  Note  1.  Many  authors  substitute  the  term  “mode”  for  “way”  when 
discussing  N-arrays.  Tucker*^  introduced  this  term  to  describe  “a  set  of  indices  by 
which  data  might  be  classified.”  The  two  terms  are  synonymous. 

Definition  A.21.  (Array  Order)  The  order  of  an  N-way  array  is  the  number  of 
dimension  into  which  the  array  extends,  i.e.,  N.  A  scalar  is  defined  as  a  zero-order  array. 
Using  Definitions  A.  13  and  A.  14  then: 
scalar  =  0-array  =  array  of  order  0, 
vector  =  1 -array  =  array  of  order  1, 
matrix  =  2-array  =  array  of  order  2, 
box  =  3-array  =  array  of  order  3, 
hyperbox  =  (iV>3)-array. 

Definition  A.22.  (Array  Unfolding)  An  N-array  (N>2)  can  be  rearranged  so  that  it  is  in 
the  form  of  a  matrix  or  a  partitioned  matrix.  When  the  different  “layers”  or  “slices’'  or 
“slabs”  of  the  N-array  are  placed  into  a  single  “plane,”  this  is  referred  to  as  the  “unfolded’ 
of  “flattened’  form  of  the  array. 

Definition  A.22.  Note  1.  The  unfolding  operation  can  be  performed  in  many  ways  and 
must  be  in  the  analyses  discussed  in  this  thesis. 

Definition  A.23.  (Outer  Product)  The  outer  product  (also  regularly  called  the  direct 
product  or  tensor  product)  of  N  vectors,  a,  b,  c,. . .  whose  elements  are  a„  bj,  Ck,..., 
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generates  an  A^-array  whose  elements  are  aibfk. ...  Outer  product  of  vectors  is  indicated  by 
®  ,e.g.,a®b<8)c®.... 

Definition  A.24.  (Outer  Product  Array)  An  outer  product  array  is  any  array  which  can 
be  expressed  as  an  outer  product. 

Definition  A.24.  Note  1.  The  outer  product  of  two  vectors,  a  and  b,  is  the  outer  product 
2-array,  M: 

Uj  ajAj  ^1^3  a^b^ 

M=a®b  =  ab'=  [A,  hj  ^3  ^4]=  ^2^1  ^2^2  ^2^3  ^2^4  (A.5) 

a^bi  a^b^  a^b^ 

Clearly,  M  is  a  rank  1  matrix  since  all  of  the  colunms  (and  rows)  are  linear  combinations 
of  the  other  columns  (rows).  A  rank  1  matrix,  i.e.,  outer  product  2-array,  is  referred  to 
as  a  dyad. 

Definition  A.25.  The  outer  product  of  three  vectors,  a  0  b  0)  c ,  generates  an  outer 
product  3-array  which  is  referred  to  as  a  triad.  The  general  case  of  an  outer  product  A^-array 
produces  an  N-&d.  The  outer  product  of  three  vectors  and  the  3-array  generated  is  displayed 
in  Figure  A.1  and  represented  algebraically  in  Equation  A6.  Unfolding  any  N-array  is  a 
straightforward  extension  of  the  three-way  model  shown  in  Figure  A.1. 


144 


a  0  b  0  c  = 


>•  b 


Slabl  Slab  2  Slab  3 

Figure  A.  1 .  Unfolding  of  outer  product  3-array.  The  outer  product  of  three  vectors  a,  b, 
and  c  with  dimensions  5, 4,  and  3,  respectively.  The  result  is  a  3-array  with  3  (=  c 
dimension),  5  by  4  slabs.  Unfolding  yields  the  5  by  12  matrix. 

Figure  A.  1,  written  in  matrix  form,  is 

0 

ci^b^c^  d^b^c^  d-^b^c^  d^Jb^c^  d^b^c^  d^b^c^ 

i^b^c^  d^b^c^  d^bi^c^  d^b^c^  d^b^c^  d^b^c^  d^Jb^c^  di^ijb^^c^  d^b^c^ 

d^b^c^  d-^b'^c^  d'lJb^C'^  d^b^c^  d^b^c^  d^j^c^  d^^b^c^  d^Jb^c^ 

x^Jb^c^  d^b^c^  d^b^c^  d^Jb^C'2^  d^b^c^  d^b^c^  d^b^c^  d^b^c^  d^b^c^ 

X^b^C^  d^b^C^  d^b^C^  d^b^C'2^  d^b^^^ 

(A.6) 

DeHnition  A.26.  (Multiple  Product)  The  multiple  product  of  N  matrices,  all  having  the 
same  number  of  columns  R,  is  the  sum  of  their  R  N-ads,  specifically, 

®b,  0c,0- 

r=l 

is  the  multiple  product  of  >^bxR’^cxR  >’"• 


Unfold  operation 
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Definition  A.27.  (A-adic  Decomposition)  For  the  N-array,  X,  the  A^-adic  decomposition 
is  the  set  of  matrices  each  having  R  columns,  whose  multiple  product  is  equal 

toX,  i.e., 

Definition  A.28.  (Array  Rank)  The  rank  of  an  A^-array,  X,  is  the  minimum  value  of  R 
(see  Definition  A.27)  for  which  an  A^-adic  decomposition  of  X  exists,  i.e.,  the  minimum 
number  of  AT-ads  which  is  equal  to  X.  The  array  rank  of  X  is  indicated  by  j?  (X) . 

Definition  A.28.  Note  1.  By  definition,  all  AAads  are  rank  1  arrays. 

Definition  A.29.  (Kronecker  Product)  The  Kronecker  product  (also  called  the  direct 
product)  of  two  matrices  is  formed  by  generating  the  product  of  each  element  of  one  matrix 
with  each  element  of  the  other  matrix  in  an  ordered  way.  The  Kronecker  product  of 
and  is  ,  defined  as: 


C=A(8)B  = 


fljjB  UjjB  •••  flj^B 

02, B  O22B  •  02*B 


[a„,B  o„2B  •••  o^Bj 

Definition  A.30.  (Positive  Definite  Matrix)  A  symmetric  matrix  is  positive  definite  if 
and  only  if  all  of  its  eigenvalues  are  greater  than  zero. 

Definition  A.31.  (Positive  Semidefinite  Matrix)  A  symmetric  matrix  is  positive 
semidefinite  if  and  only  if  all  Of  its  eigenvalues  are  greater  than  or  equal  to  zero. 
Definition  A.31.  Note  1.  Positive  semidefinite  matrices  are  also  called  nonnegative 
definite  matrices. 
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Definition  A.32.  (Gramian  Matrix  or  Gramian  Association  Matrix)  A  positive 
semidefinite  matrix. 

DeHnition  A.32.  Note  1.  There  are  four  types  of  Gramian  matrices.  These  are  formed 
from  the  data  matrix,  A,  by  taking  A 'A  or  AA'  after  any  preprocessing: 

( 1 )  Inner  product  matrix  (Data  unstandardized) 

(2)  Cosine  matrix  (Data  scaled  to  unit  variance) 

(3)  Covariance  matrix  (Data  centered  about  mean) 

(4)  Correlation  matrix  (Data  standardized) 

See  Chapter  3  for  descriptions  of  centering,  scaling,  and  standardizing. 
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APPENDIX  B.  THREE-MODE  MODELS 


Tucker  Models,  Three-Mode  Factor  Analysis,  and  TrUinear  Models 

Analysis  of  iV-arrays  was  pioneered  by  Tucker^^  who  developed  the  unfolding 
methodology  for  7V-arrays  discussed  in  Appendix  A.  His  notation  and  methodology  has 
been  the  backbone  of  multimode  analysis  ever  since.  Therefore,  a  brief  explanation  of  his 
methodology  is  warranted.  The  discussion  will  begin  with  the  2-array  notation  and  then 
progress  to  the  more  complex  3-array  notation. 

B.l  Tucker  Models 

Tucker’s  notation  involves  using  pre-subscripts  and  post-subscripts  around  the  letter 
representing  the  matrix  to  designate  the  mode  of  the  row  and  column  information, 
respectively.  Choosing  the  lower  case  italic  letters  i  and  j  to  represent  the  two  modes,  a 
two-mode  matrix  (2-array)  oriented  i  by  j  would  be  written  X^. .  The  subscript  has  several 

related  but  distinct  roles:  (1)  as  a  general  identification  of  the  mode,  (2)  as  a  subscript 
identifying  the  mode  to  which  an  element  belongs,  and  (3)  as  a  variable  identification 
symbol  for  the  elements  in  the  mode.  For  illustrative  purposes,  let Xj  be  a  matrix  of 

fluorescence  intensities  for  a  mixture  of  fluorophores  measured  at  different  excitation  and 
emission  wavelengths,  i.e.,  an  EEM. 

The  first  usage  of  this  notation  would  be  to  designate  mode  i  as  fluorescence  excitation 
wavelengths  and  mode  j  as  fluorescence  emission  wavelengths.  The  second  usage  is 
illustrated  by  assigning  identification  symbols  to  the  different  excitation  wavelengths,  i.e., 

1, ,  2, ,  3, ,  •  •  •,  AT, .  The  numeral  in  the  identification  symbol  is  the  index  value  of  the 
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element,  and  the  subscript  indicates  the  mode.  So,  if  one  measured  fluorescence  intensities 
as  a  function  of  excitation  at  10  nm  increments  from  250  nm  to  350  nm,  3,.  would 

represent  270  nm  and  as  an  index  for  any  intensity  at  270  nm.  The  number  of  elements  in 
mode  i  is  given  by  Ni  and  in  mode  j  by  Nj .  Thus,  the  dimension  of  ,  Xj  is  by  Nj . 

The  third  use  of  the  mode  symbol  is  as  a  variable  identifier  of  a  scalar  element  of  the  2- 
aiTay,e.g.,  Xy. 

Using  this  notation,  the  transpose  of  the  2-array  ,  X^. ,  is  simply  written  as  j  X,. . 

Writing  the  product  of  two  2-arrays  joins  the  post-subscript  of  the  left  2-array  with  the 
pre-subscript  of  the  right  2-array,  since  they  must  be  identical,  e.g., 

or. 

Now  consider  the  factor  analysis  of  the  EEM, X^ . 

Factor  analyzing  the  EEM,  ,  X^. ,  into  three  matrices  can  be  expressed  as; 
,X,.=.A„G„B.+N  (B.l) 

where  ,  A„ ,  is  an  N,  by  N„  orthonormal  factor  matrix  whose  columns  are  the 
eigenvectors  of  X^X, ,  ^  is  an  Nj  by  orthonormal  factor  matrix  whose  columns 
are  the  eigenvectors  of  j  X,  X^. ,  „  G„  is  a  diagonal  matrix  of  square  roots  of  eigenvalues 
associated  with  the  Gramian  inner  product  matrices  ,  X^X,  and  j  X,X^ ,  also  called  the 
“core”  matrix,  and  N  is  a  measurement  error  matrix  of  dimension  N,  by  Nj .  This  is 
equivalent  to  the  singular  value  decomposition  of  ,  X^ ,  i.e.. 
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X=USV'. 


(B.2) 


The  number  of  singular  vectors,  N„ ,  is  chosen  to  equal  the  rank  of  .  Therefore,  after 
selecting  only  the  first  N„  singular  vectors  and  values, 

.•A„G„B,=,U„S„V,-.  (B.3) 

in  njinnj  ^  ^ 

Physically  meaningful  factors  are  obtained  by  coordinate  rotation  of  the  eigenvectors  by  a 
square  transformation  (or  rotation)  matrix  T.  Letting  ,F„=,.  A„G„ ,  then, 

,F..=,r.T,.  (B.4) 

(b.5) 


where  n*  is  the  transformed  derivational  mode  and  „T„,T„  ’  =  I . 

n  n*  n 

The  extension  into  the  three-mode  case  is  straightforward;  however,  it  involves  an 
adjustment  to  the  notation  to  accommodate  higher  order  arrays.  Consider  a  3-array  with 
modes  i,j,  and  k,  X.  The  elements  in  the  array,  ,  form  a  box,  which  has  the  dimensions 

Ni ,  Nj ,  by  N^.  Now,  let  the  array  be  unfolded  so  that  it  now  is  a  matrix  with  elements 


ordered  i  by  (jk),  where  (jk)  is  read  as  j-outer  loop,  k-inner  loop.  Such  an  arrangement  is 


i  '^Uk)  ~ 


^(Jk) 

For  example. 

let 

^(jk) 

be  a 

5  by 

Nj  = 

4  by 

=3 , 3-array,  such  that 

^111 

*^112 

■^113 

^121 

•^122 

•^123 

*^131 

^132 

^133 

Xl4l 

^142 

^143 

-^211 

•^212 

■^213 

•^221 

^222 

•^223 

^231 

•^232 

•^233 

•^241 

■^242 

•^243 

^311 

•^312 

•^313 

^321 

^322 

^223 

*^331 

^332 

■^333 

^341 

*^342 

•^343 

.  (B.6) 

■^411 

•^412 

•^413 

•^421 

^422 

•^223 

•^431 

•^432 

^433 

^441 

•^442 

•^443 

-*^511 

^512 

^513 

•^521 

•^522 

^223 

•^531 

^532 

•^533 

■*541 

•^542 

*^543. 

The  notion  of  y-outer  loop,  k-inner  loop  is  easily  visualized  by  the  following  diagram: 
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AAA. 


12  3 

'■k 


1,  2,  3,  1,  2,  3, 


<X,,=2 


Based  on  the  discussion  thus  far,  note  that 

2)  The  unfolding  method  illustrated  in  Appendix  A  represents  an  i  by  {kj) 
arrangement  (considering  i,  k,  and  j  to  represent  the  a,  b,  and  c  modes,  respectively). 
Two  other  important  array  orientations  are  presented  below.  These  will  be  used  in  the 
factor  analysis  in  this  thesis: 


;Xo.)  = 


•^111 

^112 

•^113 

^211 

*^212 

•^213 

^311 

•^312 

^313 

^411 

^412 

^413 

•^511 

^512 

•^513 

•^121 

•^122 

-^123 

^221 

^222 

^223 

•^321 

^322 

•^323 

■^421 

•^422 

X423 

■’‘•521 

^522 

•^523 

*^131 

•^132 

•^133 

•^231 

•^232 

^233 

^331 

•^332 

*^333 

^431 

•^432 

^433 

^531 

^532 

•^533 

„^141 

“^142 

•^143 

^241 

•^242 

^243 

•^341 

•^342 

^343 

^441 

^442 

■^443 

^541 

^542 

*^543  „ 

(B.7) 


*  X(y)  - 


^111 

•^121 

•^131 

■^141 

^211 

•^221 

•^231 

■^241 

*^311 

•^321 

•^331 

■^341 

^411 

^421 

•^431 

•^441 

^511 

^521 

^531 

_ 1 

*^112 

•^122 

•^132 

^142 

^212 

•^222 

•^232 

^242 

^312 

•^322 

•^332 

*^342 

^412 

-^422 

^432 

♦^442 

•*-512 

^522 

•^532 

■*•542 

3i3 

^123 

^133 

•^143 

•^213 

^223 

^233 

^243 

^313 

^323 

^333 

■^343 

^413 

^423 

^433 

^443 

•*"513 

^523 

^533 

•*^543. 

(B.8) 

B.2  Three-Mode  Factor  Analysis 

Three-mode  factor  analysis  involves  decomposing  the  3-array  into  orthonormal  factor 
matrices  and  a  core  array.  In  summational  form,  the  array  is  written: 

A  g^g+f^ijk  (B-9) 
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where  a,  b,  and  c  are  the  elements  of  the  three  factor  matrices  and  g  is  an  element  of  the 
three-mode  core  matrix.  In  matrix  notation, 


(B.IO) 

where  the  notation  is  the  same  as  the  2-array  case  except  for  ,  which  is  an  unfolded 

core  3-array.  The  matrices  ,  and  are  extracted  from  the  Gramian  inner 

product  matrices, 


(B.ll) 


(B.12) 

(B.13) 


respectively,  by  solving, 

.■m,a„=:a„  (B.14) 

V,  (B.16) 

where  ,  and  v„  are  the  eigenvalue  matrices  associated  with  A„ ,  ^ B ^ ,  and  , 

respectively.  significant  eigenvectors  of  ,.M,.  are  used  to  form  ,  A^ ,  significant 
eigenvectors  of  ,.P,  are  used  to  form  ,.B„ ,  and  significant  eigenvectors  of  ,.Q,.  are  used 
to  form  C„ .  The  core  array  is  then  given  by 


(b.i?) 

where,  „A^  etc.  Since  ,A„,  ^Bp,and  are  all  column-wise  sections  of 

column- wise  orthonormal  matrices,  equation  B.17  becomes  simply 
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(B.18) 


Tucker’s  approach  is  only  exact  if  all  of  the  eigenvectors  of  the  Gramian  matrices  are  used, 
i.e.,  if  i=m,  J=p,  and  k=q.  Tucker  noted  that  “these  procedures  do  not  produce  a  least- 
squares  approximation  to  the  data.”  Kroonenberg  and  De  Leeuw®^  corrected  this  deficiency 
by  developing  an  alternating  least  squares  (ALS)  algorithm  which  they  called  TUCKALS3. 

The  TUCKALS3  solves  iteratively  for  values  of  ,A„,  and  in  turn. 


using  the  new  estimates  for  these  matrices  in  subsequent  steps. 

A  Substep: 

.  M  .  —  .  X  (y^)  (  y  B  (8)  ^  Cj;.  )  (yjt )  X  . 

(B.19) 

(B.20) 

B  Substep: 

(B.21) 

B  =  P  B  (  B  P  P  B 

J  p  J  J  p^p  J  J  J  p^ 

(B.22) 

C  Substep: 

it  Q*  ®  J  ®  p®7  ^ (y) 

(B.23) 

(B.24) 

Starting  values  of  ^  A„, ,  ^  ,  and  ^  are  obtained  from  Tucker's  method  in  Equations 

B.12-B.14.  After  evaluating  the  substeps,  the  core  array,  ,  is  estimated  using 

Equation  B.16,  completing  an  iteration.  After  each  iteration,  the  fit  is  evaluated  by 
examining  the  residuals  of  the  data  and  the  fit.  Termination  is  based  on  reaching  a 
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minimum  of  the  square  of  the  norm  of  the  residual  array,/,  i.e., 


(B.25) 

This  may  also  be  written, 
/  = 


(B.26) 


„  is  a  measure  of  the  amount  of  the  system  variance  explained  by  the  model.  The 


squares  of  the  elements  of  the  core  array,  gl^,^ ,  indicate  the  amount  of  variance  explained 

by  the  corresponding  combination  of  factors  (m,  p,  and  q)  from  each  mode. 

Transformation  of  the  three-mode  principal  factors  is  performed  using  three 
transformation  matrices,  ,  and  . 


A  =  A  T 

i^m*  i  m  m* 


(B.27) 


(B.28) 


(B.29) 


The  transformed  core  matrix  is  then 


m*^m^  (pq)^  p*^q  q 


(B.30) 


It  is  possible  to  manipulate  the  core  matrix  in  three-mode  analysis  to  create  a  more 
parsimonious  model.  This  may  be  done  by  setting  small  values  of  to  zero,  during  or 

after  three-mode  analysis,  or  by  using  chemical  knowledge  of  the  system  at  hand  to  select 
elements  of  the  core  matrix  which  are  nonzero  and  setting  the  rest  to  zero.  A  special  case 
of  the  restricted  Tucker  model  is  the  trilinear  model. 
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B.2  Trilinear  Models 


The  trilinear  model  has  core  modes  (m,  p,  and  q)  of  equal  dimension,  so 

.  and  the  core  array  is  superdiagonal,  i.e.,  the  ith  diagonal  element  of  the  ith 

slab  of  the  core  array,  is  nonzero;  and  all  other  elements  are  zero.  For  example,  the 
unfolded  core  array. 


m^(pq) 


a 

0 

0 

0 

0 

0 

0 

0 

O' 

0 

0 

0 

0 

b 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

c 

(B.31) 


with  a  ^0  and/or  b^O  and/or  c  9^  0  is  superdiagonal.  In  trilinear  analysis,  one  seeks  to 
find  the  minimal  decomposition  of  a  3-array,  i.e.,  the  smallest  value  of  needed  to  fit  the 

data  or  to  find  N^=R  (,  ) .  The  core  array  is  usually  not  utilized  explicitly  when 

performing  trilinear  analysis,  rather  the  model  is  represented  as  an  A^-adic  decomposition. 


X  =  ®(A,B,C)  +  N.  (B.32) 

Solutions  for  A,  B,  and  C  are  obtained  by  an  ALS  scheme,  specifically,  three-mode 
alternating  least  squares  (3M-ALS),  also  known  as  PARAFAC.  Setting  the  superdiagonal 
elements  of  the  superdiagonal  core  array  equal  to  one  reduces  equation  B.IO  to  equation 
B.32.  The  algorithm  proceeds  as  follows: 

A  Substep: 

(B.33) 

(B.34) 
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B  Substep: 


(B.35) 

pB,=,P(;>X,  (B.36) 

C  Substep: 

(B.37) 

(B.38) 

Iteration  is  complete  when  a  minimum  in/, 

/=||X-(®(A,B,C)f  (B.39) 


is  attained.  The  trilinear  model  has  the  very  desirable  property  that  its  decompositions  are 
rotationally  unique,  under  certain  assumptions,  vide  infra.  This  implies  that  any  minimal 
decomposition  of  the  3-array  X  will  generate  equivalent  decompositions,  i.e., 
decompositions  which  are  identical  except  for  scaling  and  permutation  ambiguities. 

To  explain  the  conditions  under  which  rotational  uniqueness  arises,  some  definitions 
are  necessary. 

DeHnition  B.l.  (Universal  A:-column  independence)  If  a  matrix  A  has  j  linearly 

independent  columns,  then  A  is  said  to  have  column  rank  j.  This  does  not  imply  that  all 
possible  sets  of  j  columns  of  A  are  linearly  independent.  However,  if  every  set  of  k 
columns  of  A  is  linearly  independent,  then  A  is  said  to  have  universal  k-column 
independence. 

Definition  B.l.  Note  1.  All  matrices  have  universal  0-column  independence.  If  the 
matrix  has  no  0  columns,  then  it  also  has  universal  1 -column  independence.  In 
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addition,  if  the  matrix  has  no  0  columns  and  has  no  identical  or  proportional  columns, 
then  it  has  universal  2-column  independence. 

DeHnition  B.2.  (^-rank)  The  maximum  value  of  k  for  which  a  matrix  has  universal 
^-column  independence  is  called  the  fc-rank. 

A  corollary  to  a  theorem  of  Kruskal,^^  which  follows,  states  the  condition  required  for 
rotational  uniqueness.  Let  Ja,  Jb,  and  Jc  be  the  k-ranks  of  the  three  /?-columns  matrices  A, 
B,  and  C,  respectively.  If 

Ja  +  A+/c  ^2/?  +  2  (B.40) 

then  0(A,B,C)  is  rotationally  unique. 


157 


APPENDIX  C.  EIGENANALYSIS-BASED  PROCEDURES 


Direct  Trilinear  Decomposition 

This  is  a  capsule  review  of  the  eigenanalysis-based  procedures  (EBP)  for  direct  trilinear 
decomposition  obtained  from  Reference  10.  References  36  and  62  contain  more  detailed 
derivations  of  the  method. 

Consider  the  simplest  trilinear  composed  of  two  bilinear  matrices,  say  EEMs  of 
mixtures,  each  containing  differing  amounts  of  the  same  fluorophores.  The  trilinear  model 
for  this  system  is 

X  =  ®(A,B,C)  (C.l) 

where  X  is  an  i  by  j  by  k  3-array  with  k  =  2  and  A,  B,  and  C  are  i,j,  and  k  by  r, 
respectively.  Also  assume  that  all  elements  in  C  are  nonzero. 

The  two  i  by  j  slices  of  X  may  be  written  as: 


X,  =AD,B' 

(C.2) 

X2=AD2B' 

(C.3) 

where  Dj  is  a  diagonal  matrix  whose  diagonal  elements  are  the  first  column  of  C,  and  Dj 
is  a  diagonal  matrix  whose  diagonal  elements  are  the  second  column  of  C.  Both  X,  and 
X2  share  common  column  and  row  spaces  whose  eigenvectors  are  U  and  V,  respectively. 
These  can,  in  turn,  be  related  to  A  and  B  by  transformation  matrices  P  and  Q,  respectively, 
i.e.. 


A  =  UP 

(C.4) 

B  =  VQ 

(C.5) 
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Since  it  is  a  simple  matter  to  compute  the  eigenvector  matrices  U  and  V,  the  task  at  hand  is 
to  compute  the  correct  transformation  matrices  P  and  Q. 


To  solve  for  P  and  Q,  let 

L,  =U'X,V  (C.6) 

and 

L2=U'X2V  (C.7) 

Combining  equations  C.2,  C.4,  and  C.6, 

L,  =U'X,V=U'AD,B'V=U'UPD,Q'V'V  =  PD,Q'  (C.8) 
and  equations  C.3,  C.5,  and  C.7, 

L2  =U'X2V  =U'AD2B'V  =U'UPD2Q'V'V  =PD2Q'  (C.9) 
Taking  the  product  of  LjLj'  yields 

LiL-'  =  PDjQ'Q -'D-‘P-'  =  PD,D-'P-'  (C  IO) 

and  the  product  of  ,  gives 

L;Lr‘  =  QDIPT* '-'D'  ’Q-'  =  QD,D-'Q"'  (C.  1 1) 

Now, 

L,L"2P  =  PD,D;’  (C.12) 


and 


L(Lr'Q  =  QD,D-'  (C.13) 

which  are  both  in  the  form  of  eigenvalue  equations  and  can  be  solved  directly  for  P  and  Q. 

The  extension  to  multiple  slices  in  the  third  mode  is  accomplished  by  creating  two 
matrices  usually  called  G  matrices,  i.e.,  G,  and  G2 .  The  G  matrices,  for  generalized 
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matrices,  are  linear  combinations  of  all  of  the  slices  in  the  3-array  such  that  all  of  the 
possible  dyads  are  included  for  the  analysis.  The  G  matrices  can  be  computed  in  a  manner 
similar  to  L,  and  Lj  and  are  used  in  their  place  in  the  calculations. 

After  the  estimates  for  A  and  B  are  determined,  C  is  determined  via  least  squares 
estimation  as  in  3M-ALS. 
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